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ABSTRACT 


This report describes the MITRE version of the TRACE-D nroeram 
now in operation. While the primary function of TRACE-D is orbit 
determination, options are also available in the program for trajec¬ 
tory prediction and observational data generation. A functional des¬ 
cription of these features is contained in this report along with a 
complete user's manual and a brief program description. 
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PREFACE 


The program described herein is a modified version of the TRACE-D 
program originally developed by the Aerospace Corporation in Los 
Angeles. Although the modifications made were many and included con¬ 
version of the program to FORTRAN IV and STRAP, restructuring of the 
program logic, and the addition of minor options and features, tne 
basic functions (i.e., orbit determination, trajectory prediction, 
and observational data generation) remain intact. 

The people involved in various aspects of the effort, in addition 
to the uthor, were Dr. K. K. Maitra, G. M. Hyder, P. E. Shifres, and 
S. Schwart 2 . A great deal of appreciation is due to K. W. Hubbard, 

C. C. Tonies, and other members of the Flight Mechanic’s Department 
of the Aerospace Corporation whose cooperation in obtaining the pro¬ 
gram and unlimited assistance throughout the task were invaluable. 

This report is intended solely as a user’s euide to the TRACE-D 
program, and therefore contains only a brief functional description 
of the analytics. For the reader seeking the mathematical details 
and arguments the ideal reference is the original document 1 issued 
by the Aerospace Corporation (see [1]). 


A Material from this Aerospace document has been freely incorporated 
into this report wherever pertinent and is not to be taken as 
original. 
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Section 1. 


BASIC INFORMATION 

1.1 Introduction 

This document is intended as a user's manual and therefore does 
not discuss the theory and mathematical details of TRACE-D except in 
a functional context. Section 1 of this report presents some program 
history, a brief summary of the modification process, and a statement 
of the program's present status. Section 2 specifies the functional 
capabilities of TRACE-D while Sections 3 and 4 describe input and out¬ 
put (patterned after the original TRACE-D document) and Section 5 
contains programming information. Wherever pertinent, throughout the 
report, material from the Aerospace documentation has been freely in¬ 
corporated . 

1.2 Program History 

The TRACE-D orbit determination program was developed by the 
Aerospace Corporation as a multipurpose, flexible computational tool 
for application to problems in orbit analysis. A major design objec¬ 
tive was that the program be able to provide answers to a wide range 
of problems relating to orbit and system design, space vehicle perfor¬ 
mance, and force model analysis. It was not and is_ not a "real time" 
program and therefore it is most effective for post-flight analysis 
and also for research and investigative purposes. The program has 
been used extensively and successfully by the staff at Aerospace in 
just such an environment. 

TRACE-D was originally written in FORTRAN II and FAP for use on 
the IBM 7094 computer in the "CHAIN" multi-coreload mode of operation. 
About thirty-percent of the program consisted of machine language 
(FAP) subroutines including the numerical integration package. The 
status of the program at the time of its release to us was reported 
to be operational for all functions. 
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1.3 Modification Process 


The primary objective of this effort was to produce a version of 
TRACE-D operable on the IBM 7030 computer here at MITRE. To accom¬ 
plish this, two separate yet simultaneous endeavors were carried out; 
one being to get the original TRACE-D program to run successfully on 
a nearby 7094 installation, and the other the actual conversion of 
routines to FORTRAN IV and STRAP and the reconstruction of the program. 
Upon completion of this phase, the next step was the checkout and 
"debugging" of the 7030 version of the program. This task was greatly 
facilitated and accelerated by running test cases for all options with 
both versions (7094 and 7030) of TRACE-D and checking all results. 

I ring accomplished the primary objective satisfactorily, work 
began on minor adjustments and additions to the program. Such fea¬ 
tures as observation residual plotting, an expanded earth geopotential 
capability, and a dynamic constraint matrix capability for site loca¬ 
tion adjustments were implemented. At the present time the program 
is operational on the 7030 computer and has produced meaningful and 
significant results for many MITRE orbit-related problems. 
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Section 2. 

MAJOR PROGRAM FUNCTIONS 
2.1 Traiectory Prediction 

Basic to all the functions of the TRACE-D program is the genera¬ 
tion of a time history (either in a forward or backward direction; 
see Section 3) of the position and velocity of a space vehicle in an 
inertial frame of reference. However, this inertial ephemeris as 
well as its associated ground track and altitude history are, in many 
instances, of interest in themselves. Thus, these computations are 
performed and printed out in the trajectory prediction mode of TRACE-D. 

The motion of a space vehicle is generated by numerical integra¬ 
tion of the appropriate differential equations of motion. Using a 
Cowell formulation (i.e,, expressing the total acceleration vector as 
three components in a cartesian system [x, y, z]) three non-linear 
second-order differential equations are derived, A predictor- 
corrector numerical integration procedure based on "eighth-order dif¬ 
ferences" is used to obtain position (x, y, z) and velocity (x, y, z) 
at any time t jased on the known position and velocity at some ini¬ 
tial time t o. 

The total acceleration vector, mentioned above, is actually a 
sum of the effects of the perturbing forces that comprise the TRACE-D 
dynamic model. The earth geopotential effects (gravity) are included 
in the form of a spherical harmonic expansion with provision for 
zonal harmonics J 2 through J 15 and all tesseral and sectorial 
terms through Ji 5 ^ 15 . Effects due to other bodies in the solar 
system (Sun, Moon, Venus, Mars, Jupiter) are computed from inverse- 
square law formulas and positions of the other bodies are obtained 
from tabulated coordinates stored on magnetic tape (planetary tape). 
Acceleration due to atmospheric drag is assumed to be directly pro¬ 
portional to the square of the velocity of the vehicle relative to 
the air. The atmospheric density is obtained from one of two different 
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model atmospheres that are incorporated into the program according to 
the user’s choice. Instantaneous changes in the inertial velocity 
vector may be applied at specific times to simulate maneuvers such as 
orbit adjust or vehicle separation. Also, an included low-thrust ac¬ 
celeration term may be. used to simulate thrust tailoff in cases in¬ 
volving large engines, long-term constant thrust, or, in some instances, 
leaking tanks or valves. 

Although all of the above effects are programmed into the TRACE-D 
equations of motion, the actual model that is employed may include 
some, all, or none of these perturbations. The user has the option 
to choose the model by setting a series of indicators. The compara¬ 
tive lasting of associated input and output quantities given in Table 
1 suggests the range of potential application. Complete instructions 
for preparation of trajectory mode input data and a sample of a typi¬ 
cal trajectory mode output listing are given in Sections 3 and 4. 

2.2 Observation Generation 

The object of the data generation function is to generate various 
forms of simulated measurements from a given definition of the ephem- 
eris of a vehicle and the location of observing stations on the earth's 
surface. These measurements may be any of twenty different types 
from up to one hundred different sites and may be generated at user- 
specified frequencies. Both "clean” and "noisy” data are available, 
the latter being generated by adding noise of specified variances and 
biases to the "clean" measurements. Upon user request, the visibility 
function (i.e., the rise and set times) of a vehicle with respect to 
one or more stations may be obtained without printout of any inter¬ 
mediate data. When data is generated, it may be written on magnetic 
tape or punched on cards and is formatted such that it is readily 
available as input to the orbit determination process. 

This option lends itself well to many applications, two of which 
are the following: studies of the effects of specific data imperfections 
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on orbit determination convergence behavior and the fitting of data 
using a different force model than that used to generate the data 
giving insight into the real-world problem of fitting live data with 
programmed (and always less than perfect) dynamic models. 

Since the calculations performed for this option encompass all 
the calculations of the trajectory prediction mode, the input quanti¬ 
ties include all those listed in Table 1 plus the additional informa¬ 
tion in Table 2, Available output quantities are listed in Table 2 
also, but a more specific list of available data types with definitions 
is included in Sections 3 and 4. 

2.3 Orbit Determination 

Orbit determination is the primary function of the TRACE-D pro¬ 
gram and encompasses both the trajectory prediction and observation 
generation modes. Stated in simplest form the orbit determination 
problem consists of extracting information from observations (possible 
20 types) of a space vehicle. The data is generally collected by a 
network of tracking stations on the surface of the earth. Nearly 
always, the information to be extracted includes the orbital elements 
but may include many other parameters of the dynamic model and of the 
tracking system. 

The process of extraction in TRACE-D takes the form of a general¬ 
ized least-squares differential correction procedure. More specifically, 
the goal is to determine values for a set of parameters (from the dy¬ 
namic and/or observational model) such that the differences between 
the actual input measurements and corresponding values computed from 
the model (usually termed observation residuals) will be minimized in 
a least-squares sense. The solution set may be comprised of not more 
than one hundred parameters, of which not more than sixty may be tra¬ 
jectory-related parameters (i.e., initial conditions, geopotential 
coefficients, etc.). The amount of data that may be fitted is re¬ 
stricted only by the nature of the problem and not by the program. 
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Table 1. 

Input/Output Quantities Associated with TRACE-D Trajectory Mode 
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Keplerian, anomalistic , 

nodal periods 
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Table 2. 

Input/Output Quantities Associated with 
TRACE-D Observation Generation Mode 


INPUT DATA 


OUTPUT DATA 


1. All input listed in Table 1. 

2. Site location coordinates 


1. Computed observations 


2. Rise-set times (per 
station) 


3. Visibility function indicator 


*3. Observation uncertain¬ 
ties 


4. Observation types 

5. Observation frequencies 

6 . Noise information 

7. Output mode (card, tape) 

8 . Refraction specifications 

9. Uncertainty matrix for 

adBARV at epoch 

* Is the propagation and effect of input 9 above. 
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In fact, TRACE-D is capable of fitting observational data from six 
vehicles in a single orbit determination procedure. It allows for 
estimation of the trajectory initial conditions and drag coefficients 
associated with five satellites in addition to the usual parameters 
in a single satellite run. 

The list of all possible input for the orbit determination mode 
is quite extensive as can be seen from Table 3. This is, of course, 
due to the comprehensive model programmed into TRACE-D for the pur¬ 
pose of wider and more varied applications. However, flexibility has 
not been sacrificed for completeness, so that the program is easily 
adaptable to either a non-real or more-real environment by simple in¬ 
put adjustment. The input and output capabilities of orbit determi¬ 
nation are listed in Table 3 and described in detail in Sections 3 
and 4. Table 4, which follows, is a complete list of all possible 
parameters that may be differentially corrected in TRACE-D. 
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Table 3. 

Input/Output Associated with 
Orbit Determination Mode 


INPUT OUTPUT 


1 . 

All input in Table 1. 

i. 

Chronological list of all 
observations by satellite 

2 . 

Parameters in solution set 

2 . 

Variational equation solu¬ 
tions 

3. 

Initial conditions for 
satellites 2-6 

3. 

Trajectory printouts 

4. 

Bounds for parameter 
correction 

4. 

Observation partial 
derivatives 

5. 

Station locations 

5. 

Residuals 

6 . 

Data errors and biases 

6 . 

Station and observation 
type RMS summary 

7. 

Data editing criterion 

7. 

Computed and predicted RMS 

8 . 

Refraction correction 
indicators 

8 . 

Corrections 

9. 

Iteration count 

9. 

Covariance matrix 

10 . 

Constraint matrix for 
parameters 

10 . 

Correlation matrix 

11 . 

Observations from all 
vehicles 

11 . 

Solution set 

12 . 

Tape specifications 

12 . 

Tape of residuals 
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Table 4. 
Parameter List* 


1. Initial position and velocity components in either spherical, 
rectangular, or classical element form for up to six indepen¬ 
dent space vehicles. 

2. Reciprocal of ballistic coefficients for up to six different 
satellites. 

3. Up to six velocity increments (KICKS) for one satellite. 

4. Amplitude and time constant for an exponentially decaying 
low thrust. 

5. Zonal harmonic coefficients J 2 through J^q and all 
te 'eral harmonics ^2,1 through Je,6« 

6 . Constant biases on all types of observations. 

7. Scale factors for range and range-rate observations. 

8 . Time biases (i.e., biases in reported times at which 
observations were made). 

9. Latitudes, longitudes and altitudes above sea level of 
observing stations. 


On any one run up to 100 of listed quantities may comprise 
the solution set with not more than 60 from types 1-5. 
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Section 3. 

Usage 

3.1 Preliminary Facts 

In order to run the TRACE-D program the user need not concern 
himself with the internal elements of the program 1 nor worry about 
handling a cumbersome deck of cards. The entire program is stored 
on tape and can be used simply bv attaching an appropriate data deck 
to the small basic running deck supplied by MITRE and submitting the 
job for execution. The basic running deck contains only tape usage 
information and the necessary instructions to trigger execution, 
therefore it is never necessary for the user to alter it in any way. 
Due to the varied tape requirements for each of the major functions, 
there are three versions of the basic running deck and the appropriate 
one must be used for each option exercised. 

The input to the TRACE-D program falls into two categories: 
one which follows the normal FORTRAN IV rules of input and the other 
which is interpreted by a special routine in the program. The former 
method is used for all observing-station data (only a fraction of 
total input) and needs no description, while the latter method is 
used for all remaining input and requires a brief explanation . 

The input routine accepts each piece of data from a data card 
in a field of eighteen columns in length with four distinct subfields 
as shown on the following page. 


1 The user who is curious or interested in the actual programming 
facts will find these concisely presented in Section 5. 
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FIELD-^ 1 

LENGTH 


10 



COLUMN NO. —> 1 2 


67 16 17 


18 


1. Prefix Field - The prefix field is one column in length and can 
occupy columns 1, 19, 37, or 55. The mode of the input that follows 
is dictated by what is punched (or not punched) in the prefix field. 
The following are the permissible prefixes. 


Prefix 

blank 

I 

D 

G 


H 


2, Location Field - The locat 
starts in columns 2, 20, 35 and 


Data-Modes 

Floating point number, number 
contains or implies a decimal 
point. 

Fixed point number, no decimal 
point. 

Hollerith information, data is 
not a number but a sequence of 
characters which can extend into 
the exponent field. 

Used in connection with integra¬ 
tion parameters, basic constants 
and the earth model. It indicates 
that card contains a symbolic name 
followed by a comment. 

Followed by 1 or 2 in column 2. 
This is simply a header card to be 
printed at various places on the 
output. The number 1 indicates 
columns 1-80 and 2 indicates 
columns 81-132 on the printout. 

field is 5 columns in length and 

The information punched in the 


location field is either a symbolic name or a subscript. A symbolic 
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name defines a location in core where the data is stored. In the case 
of a subscript, it is used with respect to the immediately preceding 
symbolic name. Therefore, whenever a symbolic name is the name of a 
table, a subscript of 1 is attached to that name and those subscripts 
greater than 1 indicate storage relative to the start of the table. 

All names and subscripts should be left-adjusted within the field. 

3. Value Field - The value field, 10 columns in length, starts in 
columns 7, 25, 43, or 61. This field together with the exponent field 
contains the numerical or hollerith data which is to be input. In the 
event the data is an integer (prefix = I) the exponent field is ignored. 
If the data is a floating point number (prefix column is blank, not 
punched) the decimal point, whenever it is not punched within the value 
field, is assumed to be after the last digit. Whenever hollerith data 
is input (prefix = D) the exponent field becomes an extension of the 
value field thus allowing for a total of twelve characters in the 
hollerith quantity. All data in the value field must be left-adjusted. 

4. Exponent field - The exponent field, two columns in length, starts 
in columns 17, 35, 53, and 71. This field allows for a range of 
floating point numbers from 10~ ig to 10“ 9 “. An exponent of - 10 or 
smaller must start with the letter J (i.e., J0 * - 10, J5 » - 15, 
etc.). As in the value field the characters in the exponent field must 
be left-adjusted. 

Although the input routine will accept four input fields per card, 
it is not necessary to use all four. One field used per card is just 
as correct since the rest of the card is ignored. 

In addition to the precise field format there are a few other 
rules of input preparation, 

1. The actual order of the data cards is almost immaterial, the only 
restriction being that all values with subscripts in the location field 
are placed in the input deck relative to the last previous symbolic 
name. 

2. Any input which does not apply to the immediate run or is to have 
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a value of zero need not be specified. 

3. In the case of two appearances of the same location symbol, the 
last one to appear defines the effective input. 

4. The basic input units are feet, degrees, and seconds. Therefore, 
whenever the units are not stated then they are a combination of the 
basic set. 

With the preceding basic description and the details on input 
preparation that follow, the user should have sufficient knowledge to 
use the TRACE-D program successfully. 

3.2 B asic Data 

The erm "basic data" is defined as the data that is common to 
all the principal TRACE-D functions. This includes the constants for 
the integration procedure, the necessary physical constants, function 
indicator, trajectory specifications (like epoch time, initial condi¬ 
tions) and force model inputs. 

3.2.1 Constants 

The constants constitute two lengthy lists of input and in¬ 
volve quite a bit of keypunching. Because of this, there is available 
a standard set, in card deck form, which may be used as is, or altered 
for specific entries. This facility reduces the user’s card prepara¬ 
tion time and eliminates possibility of error. 

In the event that the constants must be punched they must be 
formatted as shown on the following page. 
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3,2,1.1 Lines 1-19 Standard Set of Integration Constants 1 2 


LINE t 

PREFIX 

NAME OR SUBSCRIPT 

VALUE 

EXPONENT 


Cols. 

Cols. 

Cols. 

Cols. 


1,19,37,55 

2,20,38,56 

7,25,43,61 

17,35,53,71 

1 

G 

INTEG 



2 

I 

1 

1 


3 

I 

2 

2 


2 4 

I 

4 

1 


5 

I 

5 

1 


6 


11 

1.00002516 


7 


23 

1. 

JO 

8 


26 

1. 


9 


30 

1. 


10 


31 

.015625 


11 


32 

64. 


12 


33 

1. 

- 7 

13 

I 

34 

4 


14 

I 

35 

1 


15 

I 

36 

0 


16 


37 

.001 


17 


38 

0 . 


18 


40 

2820.1763 


19 

I 

41 

1 



1 For definitions of the constants consult Table 5. 

2 Notice that there are entries missing in the INTEG list, like the 
3rd one, for instance. This is because these entries are either 
not necessary for all options or are used by the program and not 
available to the user. This is true in both constants lists. 
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z. • J 

NE 

20 

21 

22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

36 

37 

38 

39 

40 

41 

42 

43 

44 

45 

46 


Lines 20 - 46 Standard Set Physical Constants 1 
PREFIX NAME OR SUBSCRIPT VALUE 


s . 

:0,38,56 

Cols. 

7,25,43,61 

C 


1 

.43752691 

2 

.55303935 

6 

.5 

8 

.43752691 

9 

1 . 

13 

3280.8399 

14 

57.2957795 

15 

20925738. 

16 

332951.3 

17 

.0122999 

18 

.814979 

19 

.107821 

20 

317.887 

21 

95.129 

22 

23454.865 

23 

3443.9336 

24 

20925738. 

25 

348762.3 

26 

32.174 

30 

348762.3 

31 

1.5 

32 

1.0471976 

33 

3.14159265 

34 

298.3 

35 

300000. 

36 

82505.922 


EXPONENT 


Cols. 

1,19,37,55 


Cols. 

17,35,53,71 


- 2 
- 2 

- 2 


ifinitions of constants consult Table 5. 
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Table 5. 

INTEG and C Constants Lists 


INTEG_ PARAMETERS OF THE TRAJECTORY INTEGRATION _ 

1 1 Formulation. 1 - Cowell (eqs. of motion), 2 - Encke 

(note: the Encke-T formulation is presently not 
available) 


2 

3 

4 


2 Differential equation subroutine. 1 - AMRK, 2 - COW 
angle in degrees for changing mean to true equinox 
1 Sun (4-9) are selectors for other-body perturbations 


5 

6 
7 


1 Moon If tape unit number at numb(18) is not zero, 

0 Venus the perturbations are included or omitted 
0 Mars according as the selector is non-zero or zero 


8 

9 

10 


0 Jupiter 
0 Saturn 

C(22) Sun (10-15) are interpolation scale factors for use 

using contents of C(22) , the number of earth-radii 


TT— 

12 

13 

1.001)02516 

C(22) 

C (22) 

Moon 

Venus 

Mars 

14 

15 

] 6 _ _ 

C (22) 
C(22) 

_ __ 0 

Jupiter 

Saturn 

Similar velocitv scale factors. (not now used, and 

thru 

21 

22 

0 

all zero). 

(used internally - Encke) 

23 

IE-10 

E Bar, COW subroutine truncation error control para¬ 
meter. 

24 


(used internally - Encke) 

25 


(used internally - Encke) 

26 

1. 

A, see COW writeup. 

27-29 


(not used) 

30 

1. 

Initial time step size. 

31 

.015625 

Minimum time step size (1/64) 

32 

64. 

Maximum time step size 

33 

IE-7 

Kepler equation convergence criterion. 

34 

4 

Ratio of Cowell step size to Runge-Kutta step size 

35 

1 

if 1, do not recompute perturbations for corrector, 
if 2, do. 

36 


Fla? for perturbation computation 
(used if INTEG(35) = 1). 

37 

.001 

Least squares convergence criterion (relative) 

38 

0 

Least squares convergence criterion (absolute) 

39 


(not used) 

40 

2820.1763 

Speed of light (earth-radii/minute) 

41 

1 

If non-zero, round position and velocity vectors 
for acceleration computation. If zero, truncate 
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Table 5 (cont.) 


INTEG PARAMETERS OF THE TRAJECTORY INTEGRATION 


42 

Reserved for later use 

C CONSTANTS 


1 

.0043752691 

Earth rotation rate (rad/min) 

2 

.0055303935 

Gm, earth gravitation constant (ER**3/min**2) 

3 


Option for speed of light correction 

4 


B=l-E, relative semi-minor axis of ellipsoid 



(computed in CSET) 

5 


B**2/A**2 = (1-E)**2 (computed in CSET) 

6 

.5 

Factor for decreasing bounds in L. S. solution 

7 


2*E-E**2 (computed in CSET) 

8 

.0043752691 

Atmosphere rotation rate (rad/min) 

9 

1 . 

Earth radius 

10 


Input N for N-sigma residual editor 

11 


Input scale factor for N-sigma residual editor 

12 


If input non-zero, go thru 1 st iteration twice 



for sums 

13 

3280.8399 

Feet/kilometer 

14 

57.2957795 

Angle conversion factor 

15 

20925738. 

A, earth radius in feet 

16 

332951.3 

Relative mass of sun 

17 

.0122999 

Moon 

18 

.814979 

Venus 

19 

.107821 

Mars 

20 

317.887 

Jupiter 

21 

95.129 

Saturn 

22 

23454.865 

Earth-radii/astronomical unit 

23 

3443.9336 

Nautical miles (6076.1155 ft)/earth radius 

24 

20925738. 

1/0 distance conversion factor 

25 

348762.3 

1/0 velocity conversion factor 

26 

32.174 

Go (used for CDA/W and thrust/W) 

27 


Input parameter difference for trajectory 



differencing 

28 


Input threshold for percent difference for 



trajectory differencing 

29 


Constant for doppler rate 

30 

348762.3 

Ft/sec per Er/rnin - 

31 

1.5- 

Factor for increasing bounds in L, S. solution 

32 

1.0471976 

Rad/min per deg/sec 

33 

3.14159265 

PI 

34 

298.3 

Reciprocal of E - ellipticity 

35 

300000. 

Critical altitude (ft) 

36 

82505.922 

Deg/day per rad/min 

37-39 


Direction cosines of body axis for look angle 
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3.2.2 


Functions to be Performed 


LINE // 


47 

D 

ITIN 

1234 



Line 47 contains the ordered list of all functions to be 
performed during a single run. Selection of functions is governed 
in accordance with the following code numbers. 

Code numbers Functions 

12 J -> Orbit Determination 

3 -> Trajectory Prediction 

4 -> Data Generation 

Up to twelve functions may be selected in any one run, i.e., the 
program will perform twelve specified functions in one continuous 
machine job. The example above specifies that an orbit determination 
be carried out first followed by a trajectory prediction and then a 
data generation. 

3.2.3 Traleetory Specifications 
3.2.3.1 Lines 48-54 Epoch 








48 


I 

YEAR 

1963 


49 


I 

MNTH 

6 


50 


I 

DAY 

15 


51 



TZNE 

0 


52 



HR 

12 


53 



MIN 

45 


54 



SEC 

15.5 



1 Please notice that the field labels are no longer specified 
(i.e., Prefix, Value, etc.) but they are, of course, implied 
by the divisions. 
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In the usual case, wherein the year, month, and day are input with the 
year positive, the X-axis is directed to the vernal equinoK. Alter¬ 
natively, if the year is input negative, the X-axis would be directed 
to the longitude of Greenwich. The hour, minute, and second entries 
refer to midnight of zone time. Greenwich Mean Time is Time Zone 0. 


3.2.3.2 Lines 55-61 Initial Conditions 


55 


I 

ICTYP 

2 


56 



IC 

126.1 


57 



2 

31.23 


58 



3 

89. 


59 



4 

14. 


60 



5 

22600114 


61 



6 

25117.3 



Line 55 indicates which of the ten IC types (1, 2, 9, 0) are 
entered in Lines 56 through 61. The alternative ICTYP entries are 
characterized as follows: 


IC Type 1 Earth-centered inertial cartesian coordinates 
(x, y, z, x, y, z in units of feet and feet 
per second) (see Section 3.1.1.1). 

IC Type 2 Spherical coordinates (a, 6, 6 A, R, v) in 
units of degrees, feet, and feet per second 
(see Section 3.1.1.2). In Line 14, negative 
r is interpreted as height above the earth’s 
surface in feet. In Line 15, if v is negative, 
cifcular velocity is computed and used. 

IC Type 3 Orbital elements (a, e, i, 0, to, t) in units of 

feet, degrees, and minutes (see Section 3.1.1.3). 

IC Type 4 Same as Item b above, with longitude X re¬ 
placing right ascension a. 

IC Type 5 No IC's input. The last trajectory point of the 
immediately preceding case is used. 

IC Type 6 No IC 1 s input. The corrected initial conditions 
from the last previous tracking run are used. 
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g. IC Type 8 Either Type 1 or 2 above, but in units of 

earth radii, minutes, and radians. Type 
number is entered at CPRAM (see Section 
5.1.2.2). 

h. IC Type 9 Same as Type 1 above, but in units of 

earth radii and earth radii per minute. 

i. IC Type 0 No IC f s input. For a tracking run, two 

R, A, E sets are used from the data to 
calculate a set of initial conditions 
(see Section 3.6.6). 

3.2.4 Force Model Inputs 

3.2.4.1 Lines 62-68 Drag and T Matrix 


62 



DRAG 

.01 


63 


I 

2 

1 


64 



4 



65 



4 



66 



6 

-8.6 


67 



7 

-55.5 


68 


I 

TMATX 

4 



Line 62 contains the drag parameter C])A/W in square feet per pound 
and Line 63 contains the atmosphere model specification (ARDC 59 or 
Lockheed/ Jacchia), The ARDC model will be used when DRAG(2) is 0 
and the Lockheed/Jacchia model when DRAG (2) is 1, If line 63 is a 1, 
Lines 64 and 65 contain quantities used in the Lockheed/Jacchia model. 

Lines 64 and 65 contain d^ and d 2 > respectively. 
These are values of certain constants in the Lockheed-Jacchia atmos¬ 
pheric density expressions. (Consult reference 1 for details). 

An entry in Line 68, TMATX, will cause the T matrix 
(the 3x3 matrix describing the dependence of drag force upon vehicle 
position) to be used in the variational equations in accordance with 
the following options: 
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a. 

TMATX = 0 

T matrix is not used. 


b. 

TMATX = 1 

Input 3p/3h is used with 
flattening. 

no earth 

c. 

TMATX = 2 

Input 3p/3h is used with 
flattening. 

earth 

d. 

TMATX = 3 

3p/3h is calculated with 
flattening. 

no earth 

e. 

TMATX = 4 

3p/3h is calculated with 
flattening. 

earth 


Also, with TMATX non-zero, Lines 66 and 67 should contain input values 
for 3p/3h. Lines 66 and 67 contain 3p/3h for altitudes between 76 
and 108 n. mi. and between 108 and 376, respectively* 

It sho d be noted that input values for 8p/3h must be used with the 
ARDC 1959 model. If a variable term is desired, the drag table 

option may be utilized. In this case, the drag parameter may be con¬ 
sidered to consist of the product of two terms, (CpA/W) x C^, wherein 
CqA/W is a constant which can be differentially corrected by the use 
of the variational equation and may be considered a function of 

Mach number below a certain altitude and as a function of altitude at 
points above that altitude. Alternatively, Cq may be considered a 
function of time. In either case, use of the drag table is necessary. 
When the drag table is not used, input of CpA/W into 
the DRAG location does not change usual TRACE operation. In this 
case, is automatically set equal to 1. If the drag table is used, 

the following additional inputs are required: 

a. DRAG(10) If DRAG(10) = 0, tables not used. 

If DRAG(10) * 1, Mach and altitude 
tables used. 

If DRAG(10) = 2, time table only used. 

b. DRAG(ll) Altitude above which altitude table is 

used and below which Mach table is used 
(needed if (DRAG(10) = 1). 

c. DRAG(12) = 0 Used by interpolation routine. 


DRAG(13) 
through (35) 


Altitude table, or time table if 
DRAG(10) = 2 (hi, C'(hi), h 2 , C^(h 2 ) 

. . • , h , C^(h ), 

11 , C^)(tj), t 2 , Gj^(t 2 ), • • • , t n , 

C D (t n)* °> °>- 
d. DRAG(36) = 0 Used by interpolation routine. 

DRAG(37) Mach table (stored as noted in Item c. 

through (59) above). 

3.2.4.2 Lines 69-75 Other Body Perturbations 


68 

I 

CTAPE 

7 



If perturbations due to other bodies in the solar system are to be 
included in the trajectory calculations, a planetary coordinate tape 
must be mounted and the logical-tape unit number 7 must be entered 
at CTAPE. 


69 


G 

C 



70 


I 

4 

1 


71 


I 

5 

1 


72 


I 

6 

0 


73 


i 

7 

0 


74 


i 

8 

0 


75 


I 

9 

0 



Lines 69-75 contain the C constants that indicate which 
bodies are being considered. The six bodies corresponding to C(4) 
through C(9) are the Sun, Moon, Venus, Mars, Jupiter, and Saturn. 
Entering a 1 in C(4) implies that the Sun is to be considered, while 
entering a 0 in C(9) implies that Saturn is ignored. Therefore, in 
the example above only the Sun and Moon are used in the other-body 
perturbation computation. 
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3.2.4.3 Lines 76-79 Exponential Thrust 


76 


1 

; 

1THRST 

1 



77 



2 


1 

78 



3 


1 

79 



4 


_1 


If an exponential thrust is to be used, the quantities T x in units 
of force/mass * ft/sec 2 , T 2 in units of min~ x , and t g and tf in 
seconds from midnight of epoch date must be input at THRST, THRST(2), 
THRST(3), and THRST(4) locations, respectively. 

3. .4.4 Lines 80-96 Gravity Perturbations 


80 

! 

1 I 

_{ 

I FLAG 

! 

81 

I 

28 

i ±i 


Lines 80 and 81 above contain the input quantity 
IFLAG(28) = 0. This is an indicator that must be set in order to 
clear storage for the earth model that is to be specified. 


82 


I 1 

10B 1 10 i 

83 


I 1 

1 2 1 6 


84 


I 

3 ! 10 

.. i 

j - 

85 


I 

! 4 I 6 . ! 

1 


Lines 82-85 define the number of terms to be used in the 
earth geopotential expansion. I0B(1) is the number of zonals and 
I0B(2) is the index of tesserals. Therefore, in the example above, 

10 zonals will be used and tesseral terms up to and including the 
6,6 term. Lines 84 and 85 specify I0B(3) and (4) which are the number 
of terms to be used in the V-matrix. The V-matrix is a measure of 
the dependence of the gravitational force upon vehicle position and is 
included as one of the many necessary calculations for the observation 
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partial derivatives (for differential correction). 

If these two values are set to zero the V-matrix calcu 
lations are ignored. For a derivation of the V-matrix and its impli 
cations Appendix B of Reference 1 is recommended. 


86 



OBJZ 



87 



2 

1,826 

-2 

88 



3 

2,485 

-6 


Lines 86-88 specify the zonal coefficients in the fol¬ 


lowing order: 


OBJZ (1) 
( 2 ) 
(3) 


not used 

Jg etc. 


89 



OBJT , 



90 



2 

1.556 

, -6 

91 



3 

2.26 

-6 

92 



4 

1.72 

-6 







93 


T 

i 

OBLT 



94 



2 

16.18 


95 



3 

35.62 


96 



4 

108.2 



Lines 89-96 specify the tesseral coefficients and the 
lambdas to be used, in the following order: 
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3.2.5 


OBJT 

^nm 

OBLT 

T nm 

(i) 

J21 

(1) 

t 21 

(2) 

CM 

CM 

•“) 

(2) 

T 22 

(3) 

J 31 

(3) 

T31 

(4) 

J 32 

(4) 

T 3 2 

(5) 

J 33 

(5) 

T 3 3 

(6) 

J4I 

(6) 

J 41 

etc. 



etc. 

Miscellaneous 

Inputs 




3.2.5.1 Lines 97-100 Instantaneous Orbit Adjusts 


97 



PKICK 

89020.31 


98 



2 

200. 


99 



3 

0 


100 



4 

62. 





















Instantaneous orbit adjusts are input at PKICK. Line 97 contains t* 
(time of first orbit adjust, OAi) in seconds from midnight of epoch; 
Line 98 contains K (magnitude of velocity change of OAi) in feet 
per second; Line 99 contains 0 (yaw angle for OAi) degrees, and 
Line 100 contains 0^ (pitch angle for 0A^) in degrees. 

The TRACE-D program will accommodate up to six orbit 
adjusts. Additional cards may be added as necessary for OA 2 
through OAg. 
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3.2.5.2 Lines 101-107 Extra Kicks 


101 


I 

NXK 

3 


102 



XKICK 

63456 


103 



2 

.2 


104 



3 

88721 


105 



4 

10 . 


106 



5 

101018 


107 



6 

.02 



Up to fifty fixed orbit adjusts (i.e., instantaneous 
changes of the in-track velocity component) may be input at XKICK. 

It should be noted that these orbit adjusts are not parameters for 
differential correction, but are applied in the equations of motion 
only and are independent of the PKICK inputs. The number of extra 
kicks (<_ 50) is input at NXK, and the table of times and AV values 
is input beginning at XKICK. The format is time, AV, time, AV, 
etc., in units of seconds from midnight of epoch day and feet per 
second, respectively. 

3.2.5.3 Line 108 Print Code 


108 

D 

PRCDE 

i 

2 

3 

4 

5 

6 

7 

8 

9 

10 

11 

12 


The print-code entry consists of two BCD words accommodating six 
character positions each. Entry of an X at any of these twelve 
positions will initiate corresponding outputs in accordance with the 
following: 

a. (1) Trajectory (trajectory only option). 

b. (2) Residuals (tracking only option). 

c. (3) Partials (tracking only option). 

d. (4) A T A after each iteration (tracking only option). 

e. (5) Variational equations (trajectory only option). 
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f . 

(6) 

g- 

(7) 

h. 

(8) 

i. 

(9) 

j • 

(10) 

k. 

(ID 

1 . 

(12) 

3.2 

.5.4 


Orbital elements (trajectory only option). 

Do not print station locations (tracking and 
data generation). 

Not used. 

Special trajectory prints 1 

Orbital elements at ascending nodes only 
(trajectory only option). 

Not used. 

Suppress all trajectory print except ascending 
nodes (trajectory only option) 2 

Line 109 Equinox Precession Corrections 


109 


DALFG 

.0002 



The rotational position of the earth with respect to the inertial 
system is characterized by the right ascension of Greenwich at mid¬ 
night on the day of epoch (oig) and is computed with respect to the 
mean equinox of epoch date. An additional factor for correcting to 
true equinox of epoch date optionally may be input at DALFG in units 
of degrees. 

3.2.5.5 Lines 110-113 Direction of Integration 

The TRACE-D program is capable of generating a time 
history of position and velocity in a forward or backward direction 
in any one run. 

The inputs which determine this direction are INTEG(30), 
INTEG(31), INTEG(32). If these inputs are negative then backward 
integration from epoch time is assumed by the program. If positive, 


1 Prints will occur at times of maximum and minimum altitude above 
the oblate earth, at times when the flight-path angle equals 90 
degrees and at special latitudes and longitudes if values are 
entered (see Section 3.3.6). 

2 If the option to write a binary trajectory tape (B7) has been 
selected, the writing of that tape is controlled by the PRTIM 
entries (see Section 3.3.1). 
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as in the standard case, forward integration from epoch is performed. 


no 


G 

INTEG 


in 


30 


-i. 

112 


31 


-.01 

113 




-10. 


The example shown above indicates to the program that 
backward integration must be performed starting at epoch with an 
initial time step of -1. minute (see Table 5 for the explanation of 
INTEG constants), a minimum time step of -.01 minutes and a maximum 
time step of -10. minutes. 

When performing any of the three basic options, these 
inputs may be set negative. In trajectory generation a time history 
back from epoch will be produced. In the observation generation mode, 
observations previous to epoch may be generated. In the orbit deter¬ 
mination mode, data previous to epoch may be fit. Generally this 
makes TRACE-D a more flexible program by making the position of epoch 
and the known initial conditions less rigid. 

It should be noted that both backward and forward inte¬ 
gration can not be performed in the same run but two runs would have 
to be made, each with a different direction of integration. 

3.2.5.6 Line 114 Primary Scratch Tape 


114 

I 

PTAPE 

4 


PTAPE is the number of the primary scratch tape that 
the program requires and 4 is its only permissible value. This card 
must be included in every data deck. 
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3•3 Trajectory-Only Data 

In addition to the basic data previously described, the following 
data is relevant to the trajectory generation. 

3.3.1 Lines 1-7 Print Time Vector Option 


1 


I 

PRTIM 



n 2 


I 

2 



rr 

O j 

u> 



3 



At i 4 



4 



5 



| 5 



At 2 6 



6 



t 2 7 



7 






1 




The above sequence of print times is for outputs selected by PRCDE 
entries (Line 108, Section 3.2.5.3). As many as nine sets of print 
intervals may exist (Line 2). In the case of the i t ^ 1 set, output is 
from t ^_i to t^ at intervals of At^, with all times in minutes 
from midnight of epoch date if PRTIM = 1 or from epoch if PRTIM = 0. 
Additional cards may be inserted if 3 ^ n 9. It should be noted 
that a normal print at epoch is automatic. 

3.3.2 Lines 8-15 Variational Equation Partial Derivatives 
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An X entered in any CPRAM, DPRAM, or KPRAM character position causes 
the corresponding variational equation to be solved. Printout of the 
partial derivatives will occur only if an X is entered at Character 
Position 5 in the PRCDE print code entry location. The ordering of 
entries in the CPRAM, DPRAM, and KPRAM character position boxes is as 
follows: 

a. CPRAM (Initial Condition Parameters) (Line 8) 

The first position specifies which one of three types 
of initial conditions is applicable, and succeeding 
positions indicate the particular parameters that are 
desired in each case. Ordering of CPRAM parameter 
entries for initial condition (IC) Types 1, 2, and 3 
as shown below. 


(IC 

Type) 


(either) 

1 

x 

y 

z 

X 

y 

z 

to 







(or) 

2 

a 

6 

6 

A 

r 

V 

to 







(or) 

3 

a 

e 

i 

a 

0) 

T 

to 







Ordering of CPRAM IC Parameter Entries 

b. DPRAM and KPRAM (Differential Equation Parameters) 

(Lines 9 through 15) 

The ordering of DPRAM and KPRAM parameter entries is 
shown below, wherein T^ and T 2 are the exponential 
thrust parameters and a^, K± 0y^, and are the 

OA number (1, 2, ..., 6 ), AV magnitude, yaw angle, 
and pitch angle, respectively. 

Sixty differential equation parameters is the maximum 
number which may be selected for any one run. 
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|DPRAM j 

1 

DRAG 

! 

J 2 I 

J 3 | 

d H | 

J5 

j d 6 

J 7 | 

Js 

J 9 

J 1 0 

d 21 

| 3 | 

^ 2 2 

JN ! 

J 32 

J 3 3 

J H 1 


1 d n 3 

J 44 

d 5 1 

d 52 

d 53 

d 5 H 

|5 i 

•J 55 

J 6I 

XTi 

! J 63 

! J64 

J 65 

1 1 .... 

! J 66 

i A 2 1 ' 

1 \ 2 2 

*31 

*32 

*33 

17 ! 

| 


i a, , 

! 4 L 

k 3 

V 

A . , 

s J 

A S2 

kT 

1 1 ..i 

^ A 54 

r 

A 55 

X 61 

X 62 

X 6 3 

! 9 ! 

HL 

A 6 4 

r. - 

1*65 

[7 

] A 66 

1- 

1 

1 

L^a 

Ti 

! t 2 

1 

1 






KPRAM la K 1 0 V 

1 ! l i y 

6 p, h i K <' \ Is 1 

a 3 * 

t-L\ 

K 3 


9p 3 

1- 

3 1 a j K | e v 

i u 1 4 y * 

i e p 4 i a5 1 Ks ! e y 5 K 

a 6 

i K » 

| e y 6 

| 9p 6 


Ordering of DPRAM and KPRAM Differential 
Equation Parameter Entries 

3.3.3 Lines 16-18 Traiectory Comparison Options 


n 1 

16 | j I IDIFF 

i 1 

» ITT 

NTAPE 

15 

18 ll I 

DTAPE 

14 

11 

-- ■ ; 1 1 

i i 

t=. = ..! 


Tape units and case indicators required for the trajectory differen¬ 
cing function are as follows: 


a. IDIFF 

b. 1DIFF 


IDIFF 


d. IDIFF 


0 A regular trajectory run is indicated. 

1 The reference trajectory will be written 
on the logical tape specified by NTAPE. 
if no entry is input at NTAPE, Logical 
Tape 15 will be used. 

2 The differences between the present and 
reference cases are computed and written 
on the logical tape specified by DTAPE. 
if no entry is input at DTAPE, Logical 
Tape 14 will be used. The difference tape 
specified by DTAPE is rewound at the 
beginning of the case. 

3 Conditions are the same as when IDIFF = 2 
except that the tape specified by DTAPE 
is not rewound. 
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e. IDIFF = 4 The tape specified by DTAPE is rewound at 

the beginning of the case and unloaded 
upon completion. 

f. IDIFF = 5 The tape specified by DTAPE is unloaded 

upon completion of the case. 

The significance of the foregoing options is that if a single-compari¬ 
son case is to be processed, IDIFF = 1 is used for the reference case 
and IDIFF = 4 for the perturbed case. If a series of perturbed cases 
are to be processed, IDIFF = 1 is used for the reference case, IDIFF 
= 2 for the first perturbed case, IDIFF = 3 for all intermediate cases, 
and IDIFF = 5 for the last perturbed case. 

When this option is used, a special running deck must be obtained 
which defines the necessary tapes specified above. 

3.3.4 Line 19 Revolution Number 


19 


REV 

7 



If an initial value other than zero is desired for the revolution 
number, it may be input at the REV location. This value must be re¬ 
initialized for each individual case. 

3.3.5 Line 20 Trajectory Tape Generation 


20 

I 

TTAPE 

3 



If TTAPE is non-zero, a binary trajectory tape will be generated on 
Logical Tape 19 in accordance with the following input options: 

a. TTAPE = 0 Tape will not be generated. 

b. TTAPE = 1 Tape will be rewound before generating 

but not unloaded after completion. This 
entry should be used for the first case 
when more than one case is involved. 

c. TTAPE = 2 Tape will not be rewound before generating 

and not unloaded after completion. This 
entry should be used for all intermediate 
cases. 
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d. TTAPE = 3 Tape will be rewound before generating 

and unloaded after completion. This 
entry should be used when only one 
case is involved. 

e. TTAPE = 4 Tape will not be rewound before generating 

but will be unloaded after completion. 

This entry should be used for the last 
case. 

This tape is not suitable for the trajectory-comparison option of the 
trajectory-only function. 

For this option also, a special running deck containing the proper 
tape definitions must be used. 

3.3.6 Lines 21-28 Latitude and/or Longitude Prints 


n, 21 


I | LATPR 

3 


22 



2 

10. 


23 



3 

15. 


24 



4 

20. 


n 2 25 


I 

L0NPR 

3 


26 



2 

200. 


27 



3 

100. 


28 



4 

180. 



Line 21 contains n l (n 1 <. 10), or the number of special lati¬ 
tudes at which trajectory prints are requested, and Lines 22 through 
24 contain the special latitudes. Additional cards may be added if 
4 n <_ 10. 

Line 25 contains ri 2 < 1^) • or ^he number of special longi- 

tudes at which trajectory prints are requested, and Lines 26 through 
28 contain the special longitudes. Additional cards may be added if 
4 <_ n 2 i 10. 

Note that an X must be entered in Character Position 9 of the 
PRCDE entry if either of the foregoing options are selected. 
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This concludes the description of trajectory related data. For 
data deck samples and running deck setups see Section 3.6. 

3.4 Observation Generation Data 

The following data is pertinent only to the data generation op¬ 
tion and is used in conjunction with the basic data (described in 
Section 3.3). 

3.4.1 Lines 1-4 Special Output Option 


1 


I FLAG 



2 

I 

6 

1 



If IFLAG(6) = 0, all generated data are printed. If IFLAG(6) * 1, 
rise, maximum elevation, and set times only are printed and the Data- 
Generation Specification Load Sheet II is not necessary except for 
listing of a card carrying TR in Columns 1 and 2 (see Section 3.4.5 
and 3.4.6). 


3 

I 

14 




If IFLAG(14) = 0, data are generated in time sequence until the 
available core space (bucket) is full. This output is then separated 
and printed by station in the same sequence as that of the input sta¬ 
tion cards. Further data are then generated until the bucket is again 
full, and the sort/print cycle is repeated. If IFLAG(14) = 1, data 
are printed as they are generated (i.e., in time sequence). 


4 

I 

ETAPE 

6 



If ETAPE is non-zero, a BCD radar observation tape will be gener¬ 
ated on the logical tape unit entered at ETAPE. The tape format will 
be the same as that of the tracking input data, including station 
locations and TF and TR cards. 
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3.4.2 


Lines 5-23 Vehicle Attitude Specifications 


y 5 



YAW 

ll. 


P 6 



PITCH 

90. 


r 7 



ROLL 

180. 



Vehicle attitude may be specified by inputting yaw, pitch, and roll 
angles in degrees in the manner shown above. These entries normally 
are introduced in conjunction with aspect angle computations. 


8 

YAW 



9 

811 

1 . 


10 

812 

120. 


11 

813 

.01 


12 

814 

.1 


13 

815 

0 


14 

816 

120. 


15 

817 

121. 


16 

818 

-.9 


17 

819 

-9.0 


18 

820 

0 


19 

821 

1420. 


20 

822 

1420.5 


21 

823 

360. 


22 

824 

0 


23 

825 

0 



The time history of vehicle attitude maneuvers may be specified by 
means of a table entered at YAW(811). The format of this table, 
which is used in connection with generation of radar aspect angles 
and may consist of up to three sets of five entries each, is itemized 
in Table 6. 
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The yaw-, pitch-, and roll-angle values that make up these input sets 
of angular rates are assumed to change at the rate given over the time 
interval defined by the start and stop times. Nominal orientation is 
zero yaw, pitch and roll which corresponds to the condition where the 
vehicle body axis is normal to the geocentric radius vector, the nose 
of the vehicle is in the in-track direction, and the top of the vehi¬ 
cle is in the direction of the extended radius vector. 

Vehicle attitude at time of epoch and for the case where the en¬ 
tries in Table 6 are all zero (i.e., when nothing is input) is assumed 
to be the attitude specified at YAW, PITCH, and ROLL. If gaps in the 
time entries of the table are present, the angles are held constant at 
the last computed values. 

Table 6. 

YAW(811) TABLE FORMAT 


Entry 

Description 

YAW(811) 

Start time in minutes from epoch for 
first set of angular rates. 

YAW(812) 

Stop time in minutes from epoch for 
first set of angular rates. 

YAW(813) 

Yaw rate in degrees per minute. 

YAW(814) 

Pitch rate in degrees per minute. 

YAW(815) 

Roll rate in degrees per minute 

YAW(816) 

Start time in minutes from epoch for 
second set of angular rates. 

YAW'(817) 

Stop time in minutes from epoch for 
second set of angular rates. 

YAW(818) 

Yaw rate in degrees per minute. 

etc. 

etc. 


3.4.3 Lines 24-38 D ata Error Inputs 


24 

B 

NOISE 377777777 

7 7 
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If NOISE is non-zero (a positive octal number), normally distributed 
random noise with standard deviation and mean value specified by input 
at RAPAR is added to the generated data. The entry at NOISE is used 
to start the random-number generator. 


25 


M 

RAPAR 

05,99 


























- 26 


D 

01,01 

H 

u 



A 

b'i 

A 

S 



27 



04,01 

.057 














28 


D 

01,02 

A 

L 



H 

B 

I 

A 

s 



29 



04,02 

200. 




The RAPAR array is used to indicate bias values for measurements 
to be generated. Line 25 merely indicates to the program that the 
array is of maximum size 05,99, Lines 26 and 27 specify that for 
station HU and the azimuth measurements the bias value to be added is 
.057 degrees. Lines 28 and 29 indicate that for station AL and the 
height measurement computed from AL the bias value to be added is 200 
feet. For additional entries the same two cards are repeated with in¬ 
creased subscripts. 
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30 i 

I 

SIGMA 

200 


31 

i 

2 i 

.... ] 

‘1 

1 

32 


3 



33 


* 

I.. 5 . 0 *.- 


34 

I |ISIG 

1 

| I 

35 ! 

1 

i 1 ! 2 

2 

1 

36 

I 1 3 

l - - - L. 

3 

- 

37 

1 1 i 4 . 

\ 

1 | 


Lines 30-37 contain the observation data weighting factors. For 
each SIGMA entry, a corresponding entry defining the sigma set and 
data type appears in ISIG Lines 34-37. The ISIG entries are of the 
form 1001 + K, where I is the observation set number and K is the 
data type. Ten sets, corresponding to I = 0, 1, 2, ..., 9, may be 
entered. This selected value of I is the same as the entry In 
Column 5 of the station location card. (See Section 3,4 y 5 .) The 
data type, K, must be one of those listed in Table 7. 

In the above example, range, azimuth and elevation data would 
contain noise with a sigma of 200 feet, ,1 degree and .1 degree res¬ 
pectively, while height would contain noise with a sigma of 5 feet. 


38 

' REFR | 

i 


: elevation is a chosen measurement to t 

)e generated, then REFR 


is the input quantity that allows for elevation refraction correction. 
The computed elevation, E, is altered to account for refraction, 
using either 

E 1 = E -t- n s i cot E 

if E ^ 0.1 radian, or 

. E + _L_ n sl " I" 6 80 

1000 12+1000E 6+1000E 

if E < 0.1 radian and n sl t 0. 
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9. The 


REFR contains the h s p term, wherein i = 0, 1, 2, 
appropriate value of i is entered on the station card in Column 6 
(see Section 3.4.5) Nominally n S o = 3.12 * 10“ 6 . Rise, set, and 
maximum elevation values are determined from the geometric E which 
represents the elevation before refraction correction is applied. Ad¬ 
ditional cards may be inserted at this location if necessary. 

Table 7. 

Data Types for ISIG Entries 


ata Type (K) 

Data Description 

Symbol 

1 

Range 

R 

2 

Azimuth 

A 

3 

Elevation 

E 

4 

Topocentric right 

a 


ascension 

T 

5 

Topocentric declination 

6 p 

6 

Topocentric hour angle 

HA 

7 

Geocentric right 



ascension 

a g 

8 

Geocentric declination 

6 g 

10 

Argument of latitude 

u 

11 

Cross plane 

V 

12 

Height 

h 

13 

X 

X 

14 

y 

i ^ 

y 

15 

z 

z 

17 

Range difference 

p 

18 

Range difference 

Q 

19 

Range rate 

R 

20 

R difference 

P 

21 

R difference 

Q 









3.4.4 Lines 39-45 Orbit Covariance Matrix 


39 



ATAS 

1. 

-4 

40 



3 

.5 

-7 

41 



6 

2. 

-2 
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10 

4 

-2 

43 



15 

2000. 


44 



21 

i. 



If observation uncertainties are to be calculated, a covariance 
matrix for the ADBARV elements at epoch must be input in lower tri¬ 
angular form at ATAS. The order of the elements is Row 1/Column 1 at 
ATAS(2), Row 2/Column 2 at ATAS(3). 

In addition, the following input is needed. 


45 

D 

CPRAM 

2XXXXXX 



The CPRAM variable indicates to the program that the covariance 
matrix has been input and in the ADBARV notation. (For a more com¬ 
plete explanation of the CPRAM variable see Section 3.3.2). 

3.4.5 Station Cards 

The locations of tracking stations which are associated with 
the observations to be generated must be input by means of station 
cards carrying appropriate descriptive information. 

Figure 1 on the next page shows the proper fields for the in¬ 
put that is described below: 

a. Station ID Columns 1 and 2 (ST). 

No two stations may be identified by the same symbol. The 
three codes that are not permissible as station ID f s are 
TS, TR, TF. 

b. Sigma Index Column 5 

The sigma index identifies the set of observation sigmas 
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to be applied to data that is generated from this 
station (see Section 3.4.3). 

c. Type of Refractivity Correction Column 6 

This indicator selects the type of refractivity 
correction which is to be applied to elevation 
readings. (This input depends on previous in¬ 
put. See Section 3.4.3). 

d. Latitude Columns 9-17 

The latitude of the station is entered in 
degrees. 

e f Longitude Columns 19-27 

The east longitude of the station is entered 
in degrees also. 

f. Altitude Columns 29-36 

The altitude of the station is entered in 
feet. 

g. P, Q Indicators Columns 38-39 and 41-42 

If P, Q, P, Q data (i.e., range difference and 
range rate difference between a master station 
and two slave stations) is to be generated these 
columns must contain the station ID 1 s of the 
slave stations. Also, it is necessary to enter 
a station card for each slave station but col¬ 
umns 38-39 and 41-42 are to be left blank. 

The station cards must be preceded by two cards with END in 

columns 1-3, and must be followed by a card with TS in columns 1-2. 

3.4.6 Data Limit Specification Cards 

a. Station ID (ST) Columns 1-2 

b. Data Rate Columns 9-16 

Time interval in seconds at which data for a given 
station are to be generated and testing interval 
for rise set option. 

c. Minimum Elevation Columns 18-23 

The minimum elevation at which the vehicle is 

visible to this station. 
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d. Maximum Elevation Columns 25-30 

The maximum elevation at which the vehicle is visible 
to this station. Zero value implies a 90 degree limit. 

e. Maximum Range Columns 32-40 

The maximum range in nautical miles to which the 
vehicle is visible. Zero value causes this test 
to be ignored. 

f. Start Time Columns 51-58 

Start time from midnight of start date in days, 
hours, and minutes. Zero value implies that 
epoch is start time. 

g. Stop Time Columns 60-67 

Stop time from midnight of start date in days, 
hours, and minutes. 

The set of cards must follow the TS card in the deck and them¬ 
selves be followed by a card with TR in the first two columns. 

3.4.7 Data Type Cards 

The following set of cards (which are the last in the card 
deck) specify the exact type of measurements to be generated from each 
station. 

Figure 3 is a card image. 

a. Station ID Columns 1-2 

b. Observation Types Columns 7-33 

An X in the appropriate column will produce specific measure¬ 
ments according to Table 8. 
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NOTES 

ST is the station identification code. 

Station codes must be entered in Columns 38-39 and 

4i-42 if P, Q or P, Q data are input. 

An asterisk (*) following a unit value indicates 
that a decimal point is necessary for that 
quantity. 
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Identification 
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NOTTS 

ST is the identifica¬ 
tion code. 

If P, P , Q, Q are 
selected the speci¬ 

fication of an ad¬ 
ditional station- 
identification code 

is necessary. 
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Table 8. 

Output Quantities Corresponding to 
Columns 7 through 33 of Data-Generation Type Card 


Column 

Output Quantity 

Unit 

7* 

Range 

n mi 

(ft on ETAPE) 

8* 

Azimuth 

deg 

9* 

Elevation 

deg 

10* 

Range rate 

ft/sec 

11-14* 

P, Q, P, Q 

ft/sec, ft 

15 

Azimuth rate 

deg/min 

16 

Elevation rate 

deg/min 

17 

Range acceleration 

ft/sec 2 

18 

Mutual visibility (Output will be 
a list of numbers of stations 
visible at output time. Stations 
numbered in order of input on 
station cards. Number of stations, 

8 maximum). 


19 

Geodetic latitude of vehicle 

deg 

20 

Longitude of vehicle 

deg 

21 

Surface range from station to 
subvehicle point 

n mi 

22* 

Altitude of vehicle 

n mi 

(ft on ETAPE) 

23 

Doppler rate 


24 

Look angle (Angle between a vehicle 
axis and the station/vehicle line 
of sight. The direction cosines of 
the vehicle axis in the basic iner¬ 
tial system must be entered in C(37), 
C(38), and C(39). These quantities 

deg 


* These quantities are output on ETAPE 
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Table 8 (cont.) 


Column 

Output Quantity 

Unit 


may be input as constant or the user 
may provide a subroutine (FANG) to 
compute the direction cosines at each 
output point). 


25 

Observation uncertainties* (If in- 

Same units 


verse A T A matrix for initial condi- 

as obser¬ 


tions is input at ATAS, the ADBARV 

vations 


elements are selected as parameters, 

and if an X is entered in Column 25, 

the [A t A] - 1 is updated to observation 

times and the standard deviations in 

• • • 

the quantities R, A, E, R, A, E are 
derived and printed. The uncertain¬ 
ties are only those due to the un¬ 
certainty in the ephemeris which is 
implied by the given [A^A]" 1 for the 
epoch conditions). 


26 

Angle kappa (K). (Angle between 
station line-of-sight and geocentric 
radius vectors). 

deg 

27 

Aspect angles. (Angle 1 ($>) is de¬ 
fined as the angle between the vehicle 
yaw axis and projection of the station 
line-of-sight vector in the roll plane. 
Angle 2 (0) is defined as the angle be¬ 
tween the vehicle roll axis and the 
line-of-sight vector to the station. 

deg 

28 

Signal attenuation = -40 logio R > where 

R is slant range in feet 

db 

29* 

x, y, z (Same rectangular earth-fixed 

n mi 


(X through Greenwich) geocentric quan¬ 
tities accepted as Type 5 observations 
for orbit determination). 

(ft on ETAPE) 


* These quantities are output on ETAPE 


52 









Table 8 (cont.) 


Column 

Output Quantity 

Unit 

30* 

Topocentric right ascension and 
declination 

deg 

31* 

Geocentric right ascension and 
declination 

deg 

32* 

Topocentric hour angle 

deg 

33* 

Vehicle-centered argument of 
latitude and cross-plan angle 

deg 


* These quantities are output on ETAPE 
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3.5 


Orbit Determination Data 


This section gives a description of all the pertinent data for 
a differential correction run to be used in conjunction with the pre¬ 
viously described basic data. 

3.5.1 Lines 1-14 Multiple Satellite Option 

As was mentioned in Section 2, TRACE-D will accept measurements 
for up to six satellites on any one run. The method of specifying the 
input for these vehicles is shown below. 


Year 1 


I 

SAT 2 

1964 


Month 2 


I 

2 

2 


Day 3 


I 

3 

10 


Hour 4 



4 

3. 


Minute 5 



5 

30. 


Second 6 



6 

52. 


ICTYP 7 


I 

7 

2 


IC 8 



8 

352. 


2 9 



9 

10. 


3 10 



10 

90.05 


4 11 



11 

165. 


5 12 



12 

22580632. 


6 13 



13 

25205.3 


DRAG 14 



14 

.015 



If observations for a second satellite are to be input, Lines 1 through 
14 are used for entry of epoch, initial conditions, and drag coeffi¬ 
cient for Satellite 2. Lines 1 through 3 contain the year, month, and 
day, and Lines 4 through 6 contain the hour, minutes. and seconds, 
Greenwich time. Line 7 indicates the type of initial conditions that 
may be entered in Lines 8 through 13 (Type 1, 2, 3, 4, 8, or 9). Line 
14 contains the drag coefficient (CdA/W). The input initial conditions 
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(IC) for satellites 3 to 6 are the same as those noted above except 
that the symbol SAT2 is replaced by SAT3, SAT4, SAT5, and SAT6, as 
appropriate. 

3.5.2 Lines 15-31 Parameter Specification 

The method of selecting the parameter set for the differential 
correction process is described below. 

3.5.2.1 Differential Equation Parameters 


15 


D 

CPRAM 2 

X 

X 

X 

X 

X 

X 








16 


D 

CPRAM x 














17 


D 

3 














18 


D 

5 














19 


D 

7 














20 


D 

9 














21 


D 

PSAT 2 














22 


D 

PSAT3 














23 


D 

PSAT4 














24 


D 

PS AT 5 














25 


D 

PSAT 6 














26 


D 

KPRAM 














27 


D 

3 















An X entered in any CPRAM, DPRAM, PSAT, or KPRAM character position 
causes the corresponding parameter to be differentially corrected. 

The ordering of entries in the CPRAM, DPRAM, PSAT, and KPRAM charac¬ 


ter position boxes is as follows: 

a. CPRAM (Initial Condition Parameters for Satellite 1) 

(Line 15) 

The first position specifies which one of these types 
of initial conditions is applicable, and succeeding 
positions indicate the particular parameters that are 
desired in each case. Ordering of parameter entries 
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for IC Types 1, 2, and 3 is as shown in Section 
3.3.2. 


b. DPRAM (Differential Equation Parameters) 

(Lines 16 through 20) 

The ordering of DPRAM differential equation para¬ 
meter entries is as shown in Section 3.3.2. How¬ 
ever, in this application the T^ and T 2 para¬ 
meters are associated with Satellite 1 only. 

c. PSAT2 through PSAT6 (Initial Condition and Drag 

Coefficient Parameters for Satellites 2-6) 

(Lines 21 through 25) 

The ordering of PSAT IC and drag parameter entries 
for Satellites 2 through 6 is shown below. 


' 

PSAT 2 

“2 

& 2 

62 

a 2 

^2 

V 2 

t0 2 

drag 2 

PSAT 3 

“3 


e 3 

A 3 

r 3 

v 3 

c o 3 

DRAG 3 

PSAT4 

“h 

6 4 

34 

a 4 

*4 

v 4 

C O4 

DRAG 4 

PSAT 5 

“5 

«5 

3s 

As 

r 5 

V5 

c °5 

DRAG 5 

PSAT 6 

“6 

«6 

36 

a 6 

r 6 

v 6 

C °6 

DRAG 6 


It Is important to note that only IC Type 2 may be specified for 

Satellites 2 through 6. 

d. KPRAM (Orbit Parameters for Satellite 1) 

(Lines 26 and 27) 

Ordering of KPRAM orbit adjust parameters is as 
previously shown in Section 3.3.2. It should be 
noted that the orbit adjusts are for Satellite 1 
only. 

Sixty trajectory parameters is the maximum number 
which may be selected for simultaneous solution. 
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3.5.2.2 Radar Parameters 


28 

M 

RAPAR 

05, 99 

















29 

D 

01, 01 

B 

N 



L 

A 

T 






30 


03, 01 

.01 


31 


04, 01 








32 

D 

01, 02 

C 

K 

2 

4 

R 

B 

I 

A 

s 




33 


03, 02 

500. 


34 


04, 02 

750. 







35 

D 

01, 03 

H 

u 



T 

B 

I 

A 

S 




36 


03, 03 

.001 


37 


04, 03 

.0035 



Line 28 indicates that the data listed subsequent to RAPAR on the 
load sheet will be input into a 5 * 99 matrix array which has been 
preset to zero. The columns of this array correspond to parameters, 
and the rows correspond to parameter identification (Positions 1 and 
2), bounds (Position 3), and bias estimates (Position 4) respectively. 
Row 5 is currently not used. 

Lines 29 through 31 specify the first applicable radar 
parameter. Line 29 contains the station name (Positions 1 and 2), 
pass identification (Positions 3 and 4), and parameter name (Positions 
5 through 10). Lines 30 and 31 contain the bound and initial value, 
respectively, in feet, degrees, and minutes except for R, P, and Q, 
which are in feet per second. If the parameter is station latitude, 
longitude, or altitude, the initial value is taken from the station 
location card. 
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If the pass identification character position is left 
blank, all data with the indicated station name will be used to cor¬ 
rect the parameter. If the pass identification is not omitted, only 
data that are identified by both the indicated station name and indi¬ 
cated pass identification will be used to correct the parameter. If 
the radar parameter specified is station latitude, longitude, or alti¬ 
tude, the pass identification is ignored and all data with the station 
name are used for the parameter correction. 

Lines 32 through 34 and Lines 35 through 37 specify the 
second and third radar parameters respectively. Station names, pass 
identifications, parameter names, bounds, and initial values are 
treateJ. in the same manner as the inputs in Lines 29 through 31 des¬ 
cribed above. Additional cards may be added in cases where more than 
three radar parameters are involved. Available radar parameters are 
listed in Table 9, 

Note that the total number of parameters which may be se¬ 
lected for simultaneous solution must be less than one hundred. 
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Table 9. 


PARAMETER 

SYMBOL 

Station latitude 

LAT 

Station longitude 

LONG 

Station altitude 

ALT 

Time bias 

TBIAS 

Range bias 

RBIAS 

Azimuth bias 

ABIAS 

Elevation bias 

EBIAS 

Topocentric right ascension bias 

RTBIAS 

Topocentric declination bias 

DTBIAS 

Topocentric hour angle bias 

HABIAS 

Geocentric right ascension bias 

RGBIAS 

Geocentric declination bias 

DGBIAS 

Argument of latitude (u) bias 

UBIAS 

Cross plane (v) bias 

VBIAS 

Height bias 

HBIAS 

x bias 

XBIAS 

y bias 

YBIAS 

z bias 

ZBIAS 

P bias 

PBIAS 

Q bias 

QBIAS 

Range-rate bias 

RDBIAS 

P bias 

PDBIAS 

Q bias 

QDBIAS 

Range scale factor (Kr) 

KR 

Range-rate scale factor (Kq) 

KD 
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3.5.3 Lines 38-50 Parameter Bounds 


38 



BNDS 

.5 


39 



2 

.5 


40 



3 

.1 


41 



4 

.5 


42 



5 

1000. 


43 



6 

5 


44' 



7 

.05 


45 



8 



46 



9 



47 



10 



48 



11 



49 



12 



50 



13 




A bound must be entered for each parameter selected. These 
bound entries must be in the same sequence as the parameters. For 
each iteration of the differential correction process, the change in 
each parameter is less in absolute value than the corresponding bound 
if that bound is positive, zero if the corresponding bound is zero, 
or unrestricted if the corresponding bound is negative. 
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3.5.4 Lines 51-62 Observation Sigmas 


51 



SIGMA 

100. 


52 



2 

.5 


53 



3 

.5 


54 



4 

1000. 


55 



5 

200. 


56 



6 

300. 


57 


I 

ISIG 

1 


58 


I 

2 

2 


59 


I 

3 

3 


60 


I 

4 

113 


61 


I 

5 

114 


62 


I 

6 

115 



The observation sigmas for each observation type must be input 
in the same manner as was previously well described in Section 3.4.3. 

3.5.5 Line 63 Maximum Iterations 


63 

I 

MAXIT 

4 



The symbol MAXIT must be input with an integer value to specify 
a stopping point to the run. If the differential correction process 
has not converged at the end of MAXIT iterations, the process will be 
terminated. 

3.5.6 Lines 64-66 Residual Editing 


64 


CEDIT 

3 


65 


2 

.9 


66 


3 

1 
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If CEDIT is zero, no editing is done. If CEDIT Is non-zero, 
residuals will be edited In accordance with the following: 

a. CEDIT ^ 0 Residuals greater than (input 

sigma * |CEDIT|) will be discarded. 

b^ CEDIT - 0 Residuals greater than (statistical 

sigma from the previous iteration x CEDIT) 
will be discarded. No editing is done on 
Iteration 1. Sigmas are computed for the 
first five data types encountered for each 
station. 

Line 65 represents a scale factor such that if CEDIT(2) is non¬ 
zero, CEDIT is replaced by (CEDIT * CEDIT(2)) at the end of each itera¬ 
tion. It CEDIT(2) is zero, CEDIT is not modified. 

Line 66 Is a special option wherein if CEDIT(3) is non-zero, 
Iteration 1 will be repeated with editing performed with the sigmas 
computed during the first pass through Iteration 1. This will allow 
editing to be done on all iterations with computed sigmas. 

Sigmas (rms) for the first five data types for each station are 
computed and printed at the end of each iteration regardless of the 
residual-editing option selected. 

3.5.7 Data Correction Factors 

3.5.7.1 Lines 67 and 68 Elevation Angle Refraction 


67 

1 

REFR 



68 


2 

4.0 

-6 


A table of refraction indices , which may contain 
up to ten values, may be input starting at REFR. The entry used to 
compute refraction corrections for radar elevation observations is de¬ 
termined by the type number contained in Column 6 of the corresponding 
station location card. A zero in Column 6 causes the entry at REFR to 
be used, a 1 in Column 6 causes the entry at REFER t 1 to be used, etc. 

If the table contains no entries, the value 
312.0 * 10 ” 6 , which is built in at location REFR, will be used to 
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compute refraction corrections for all data whose station location 
cards contain zero in Column 6. All other positions of the table are 
assembled as zeros. 

3.5.7.2 Range Refraction Correction 


69 


RREFC 

1 



Refraction corrections to all range observations will 
be computed and applied if RREFC is non-zero. 

3.5.7.3 Line 70 Propagation Time Correction 


70 


SLT 

-2820.1763 



The velocity to be used in calculating the observation 
time correction due to propagation time is entered at SLT in units of 

earth radii per minute. In the absence of an entry, no correction will 

be applied. If an entry is present the correction will be applied to 

times associated with R, A, E, R, h, P, Q, P and Q data only. 

3.5.8 Lines 71, 72 Data Tape Options 


71 


I 

IBCDI 

6 


72 


I 

IBINI 

1 



If the radar observation and station location information is to 
be input via a BCD tape other than the normal FORTRAN system input tape, 
the tape number must be specified at IBCDI as 11. If a binary tape 
containing compacted radar data produced by a previous run is to be 
input, IBTNI must be specified. 
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3*5.9 Lines 73, 74 True Equinox Correction 


73 



DALFG 

-.004 


74 



2 




An additive factor may be applied to the computed right ascen¬ 
sion of Greenwich at midnight of epoch day by entering the appropriate 
value in units of degrees at DALFG for Vehicle 1, at DAFLG(2) for 
Vehicle 2, etc. This entry usually is used to correct from mean to 
true equinox reference coordinates. 

3.5.10 Lines 75-78 Proximity Indicator 


4* 75 



AN0M1 


A 76 



2 


h 77 



3 


R 0 78 



4 



During an orbit determination run, an indicator may be obtained 
whenever the trajectory passes within a given distance (range) of a 
point on the surface of the earth by input of geodetic latitude (deg), 
east longitude (deg), and altitude (n mi) of the point at AN0M1 and 
the succeeding two positions and of the testing distance (the range 
from the point to the vehicle) at AM0N1(4). Testing and printing is 
done for up to three such sets of input at AM0N1, AN0M2, and AN0M3. 

3.5.11 Lines 79 on Parameter Constraint Options 

There are two different parameter constraint options available 
in TRACE-D. Both are described below in detail. 

3.5.11.1 General Constraint Matrix 

Option 1 allows the user to apply linear constraints 
to any of the variables in the parameter set, and to specify these 
constraints by means of an input matrix set up as shown below. 
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79 


I 

KNST 4 


80 


I 

BLIST 1 


81 



i 


bn 82 



i 


83 


I 

2 


84 


I 

2 


85 



1 


86 


I 

3 

L 


87 


I 

3 


88 



1 

i 

89 


T 

A 

! - i 

\ 

i 

90 


I 

4 

1 

b 44 91 


i 

i 

1 

[ 1 

92 


I 

5 j 

_ 1 _ 

93 


I 

1 


b 5 i 94 



-1 


95 


I 

6 


96 


T 

2 


b 6 ? 97 



.5 

i 


98 



5 


99 


I 

5 


c 51 100 



6 


101 

I 

7 


102 


I 

5 


103 



_i_^ 



The clearest explanation of how this constraint matrix 
works is an example. Therefore, in the case shown above it is as¬ 
sumed that n parameters are to be solved for where 
(pi, p 2 , p n ) = P. The ordering of the Pi corresponds to the 


65 























































order of the X f s for the CPRAM, DPRAM, and KPRAM and the RAPAR ar¬ 
rays, respectively. Further assuming that these parameters are to be 
subjected to m linear constraints, which, for example for n - 6 
and m = 2 might be Pi + P 5 = 6 , P 2 - 2p^ « 0, KNST above would be 
equal to (n - m) = 4, or the number of effective (unconstrained) 
parameters. 


The BLIST constraint matrix is represented by the factor 
B in the expression p ** Bp + c, where the p are the effective 
parameters. In connection with the foregoing example, this expression 
would assume the form 







— 


f— — 


- — 

Pi 


1 

0 

0 

0 


Pi 


0 

P2 


0 

1 

0 

0 


p 2 


0 

P 3 


0 

0 

1 

0 

X 

P 3 

+ 

0 

P4 


0 

0 

0 

1 


P4 


0 

P 5 


- 1 

0 

0 

0 




6 

P6 


0 

.5 

0 

0 




0 







v 

— — 


_ - 


V_ 



J 

V 

— 

'—' 

P 


B 




p 

c 


Pertinent constraints are applied by input of the augmented (n + 1) 
by (n - m + 1 ) matrix 



in the BLIST array. If a^j is an element of this augmented matrix, 
then i, j, and a^j would be input. The inputs for the example 
given above would be as indicated in Lines 79-103. KNST (line 70) is 
input as 4 and then the BLIST matrix follows. Lines 80, 81, and 82 
specify that BLIST (1, 1) is a 1 and lines 83, 84, and 85 that BLIST 
( 2 , 2 ) is a 1 , etc. Elements of the matrix not input are assumed to 
be zero. In this way any subset of the set of parameters may be 
linearly constrained and the constraint matrix may be inDut to the 
program in the exact same manner. 
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3.5.11.2 Special Application to Station Location 

The other option available is a specific case of the 
general constraint matrix which is applicable in the following type of 
problem. Consider a differential correction run with data from a tri¬ 
static system and a set of parameters including the latitudes, longi¬ 
tudes and altitudes of the three stations. In such a case, by choosing 
this constraint option^the correction of the station coordinates can 
be performed such that the baseline distances are held constant. 

Mathematically, this implies that three constraining 
equations 1 are applied to the nine station locations; more specifi¬ 
cally they are of the following form. 


A<J >\ - ai A 0 i + a2 AHi 4 - a.% A02 + a4 AH2 + a5 A3 4 * a$ AH3 


A 4>2 = bi A 0 } + b2 AHi 4 - b3 A02 + b4 AH2 4 * b5 A 4 >3 4 * b$ AH3 


A03 = ci A@i + C2 AHi 4 - C3 A 02 4 C4 AH2 4 * C5 A$3 4 * c$ AH3 

where 


Oi, *1, Hi) 

(02, 4>2> H 2 ) 
(©3) 4>3> H 3) 


are the three sets of latitude. 
Longitude, and altitude for 
Stations 1, 2, and 3. 


and a^, b^, c^ i = 1, ..., 6 are the 18 coefficients of the 

three equations. 


With the equations in this particular form, it seems 
that the general constraint option could be applied easily. However, 
there are two reasons for the separate treatment of this application 
that are not immediately apparent. First, the constants a^, bj_, c^ 
i = 1, ..., 6 are functions of the current station location 2 and 
are not easily and accurately calculable by the user. Second, since 


1 

The detailed derivation of this particular form of the constraining 
equations is not given here but may be found in Reference 3. 

2 

The exact expressions for the coefficients can be found in Reference 3. 
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the constants are functions of the current station locations, they 
must be reevaluated every time the station locations change (i.e., 
after every iteration of the correction process the constants must be 
reset since the station locations have been adjusted). In other words, 
the type of constraint matrix needed here is a dynamic one that is 
adaptable to change throughout the estimation procedure. It is now 
obvious that the general constraint matrix of Section 3.5.11.1 is too 
rigid for this purpose. 

Thus, added to the TRACE-D program is the capability to 
compute these necessary constants from input information, insert them 
into the appropriate location of the constraint matrix, recompute the 
const <~s as a function of the new station locations after each itera¬ 
tion and adjust the constraint matrix accordingly. Once the user has 
set up the appropriate input the rest is performed by the program. 

The input necessary is best described with an example. 
Consider a differential correction run with data from three stations 
and a parameter set consisting of the ADBARV elements, drag, and sta¬ 
tion locations. To employ the fixed baseline distance constraints on 
the station locations consider B in the following expression. 

Ap' = B AP + c 

There are 16 parameters to be estimated but there are 3 constraint 
equations, therefore there are 13 effective or unconstrainted para¬ 
meters. B is then of dimension 16 x 13 and looks like this: 


68 


Aa 


r 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0~ 


Aa 


[o“ 

A6 


0 

1 

0 

0 

0 

0 

0 

0 

0 

n 

0 

0 

0 


A6 


0 

ab 


0 

0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 


AB 


0 

aa 


0 

0 

0 

1 

0 

0 

0 

0 

0 

0 

0 

0 

0 


AA 


0 

AR 


0 

0 

0 

0 

1 

0 

0 

0 

0 

0 

0 

0 

0 


AR 


0 

AV 


0 

0 

0 

0 

0 

1 

0 

0 

0 

0 

0 

0 

0 


AV 


0 

A drag 


0 

0 

0 

0 

0 

0 

1 

0 

0 

0 

0 

0 

0 


A drag 


0 

A lat i 

_ 

0 

0 

0 

0 

0 

0 

0 

a l 

a 2 

a 3 

a 4 

a s 

a 6 


A longi 

+ 

0 

A longi 


0 

0 

0 

0 

0 

0 

0 

1 

0 

0 

0 

0 

0 


A alt i 


0 

A alt i 


0 

0 

0 

0 

0 

0 

0 

0 

1 

0 

0 

0 

0 


A long2 


0 

A lat2 


0 

0 

0 

0 

0 

0 

0 

bl 

b 2 

b 3 

b 4 

b 5 

b 6 


A alt2 


0 

A long2 


0 

0 

0 

0 

0 

0 

0 

0 

0 

1 

0 

0 

0 


A lat 3 


0 

A alt2 


0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

1 

0 

0 


A alt 3 


0 

A lat3 


0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

1 

0 


0 

A long 3 


0 

0 

0 

0 

0 

0 

0 

Cl 

c 2 

c 3 

c 4 

C 5 

c 6 


0 

A alt3 
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0 

quantities 

a i > 

t>i, 

c i> 

i 


1, 

.. 

* > 

6 , 

are 

* 

of course, 

the 


ficients in the constraint equations which are unknown to the user. 

This B matrix is input in exactly the same manner as the 
general constraint matrix except that each element that is a coeffi¬ 
cient (either ai, b^, c^) must be input with a value of - 10,000. 

This is merely to indicate to the program the position in the matrix 
of these terms. 

Lines 104 on show the appropriate format for the necessary 
input of this case. 
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Lines 104 and 105 specify the quantity NUMB(44) which must be input 
as a - 1 to choose this constraint option. Line 106 is the number of 
unconstrained parameters in the parameter set. Lines 107 on are the 
non-zero elements of the B matrix entered in order. 

3.5.12 Station Identification Cards 

The locations of tracking stations or of points on the surface 
of the earth which are associated with observations must be input by 
means of station cards carrying appropriate information indicating 
the manner in which the corresponding observations are to be processed. 
These cards were previously discussed in Section 3.4.5 and shown in 
Figure 1, but the format is included again here for completeness. 

a. Columns 1 and 2 (ST) : Station identification symbol. 

No two stations may be identified by the same symbol 
or any one station by the symbol TS. 

b. Column 5: Sigma index identifying observation-sigma 
set to be applied to data from corresponding station. 

The sets of sigmas are input with the base data. 

c. Column 6: Type of refractivity correction to be used 
for elevation readings from this station. Refracti- 
vities are numbered in their input order within the 
base data (see Line 67, REFR). 
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d. Columns 9 through 17 : North latitude of station in 
degrees. 

e. Columns 1 9 through 27 : East longitude of station in 
degrees. 

f. Columns 29 through 36 : Altitude of station in feet. 

g. Columns 38/39 and 41/42 : If a station reports P, 0, 
or £, Q data. Columns 38/39 and 41/42 contain the 
two letter symbols for the associated station(s) of 
the tracking configuration. Each such associated 
station must be represented by a separate station 
card, but it is not necessary for Columns 38/39 
and 41/42 to be filled out on the latter. 

These cards must be preceded by a card with END in Columns 1-3 and 

followed by a card with TS in Columns 1 and 2. 

3.5.13 Observation Data Cards 

The observations to be processed must be the last set of cards 
in the deck and must follow the station location cards. 

Since the number of observations that may be input to TRACE-D 
is unlimited, there is a flocking procedure that must be followed. A 
flock of data consists of 200 or less observation cards for a single 
satellite. There is no limit to the number of flocks that can be 
stacked; however, each flock must be followed by a card with TF in 
Columns 1 and 2 and the times associated with the observations of a 
flock can never be earlier than the latest time of the previous flock. 
This means that the data is not time sorted within a flock but the 
flocks are sorted according to earliest and latest times therein 
contained. 

The contents of an observation card is formatted in the fol¬ 
lowing way, as shown in Figure 4. 

a. Columns 1 and 2 (ST): Station Identification, symbol 
which must correspond to a station location card. 

b. Columns 3 and 4: Pass identification (optional). 

c. Columns 5 through 21: month (columns 5 and 6), day 
(columns 7 and 8), hours (columns 9 and 10) 
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minutes (columns 11 and 12), and seconds (columns 
13-21), of the corresponding observations relative 
to Greenwich Mean Time. 

d. Column 22: Observation set number. 

e. Columns 23-76: Observation 1 (columns 23-40), 

Observation 2 (columns 41-58), Observation 3 
(columns 59-76). 

f. Column 77: .Card number indicating observations, 
variances, or covariances. 

The observation set (d) number and the card number (f) are the two in¬ 
dicators that define exactly what type of observations are contained 
on each card. Table 10 gives a complete listing of all types according 
to observation set number and card number. 

After the last observation card in the deck, a TF card is not 
necessary, however a card with TR in Columns 1 and 2 plus an END card 
must be attached as the last two cards of the deck. If data from more 
than one satellite are to be used, the Satellite 1 data are set up in 
the same manner except that the TR card must be replaced by a TT card. 
The data for Satellite 2 are then similarly arranged. If Satellite 2 
is the last vehicle, corresponding data are followed by a TR card and 
an END card: if it is not the last satellite, data are followed by a 
TT card. Data for Satellites 3 through 6, as applicable, are added in 
the same manner. A TR card rather than a TT card must follow the data 
corresponding to the last satellite, and an END card must follow the 
TR card. See Section 3.6 for a picture of the correct setup. 

3.6 Input Deck Arrangements 

Since the basic running deck to be used for each function to 
be performed can be obtained, the user need not concern him¬ 
self with its contents. Therefore in each of the Figures the separate 
parts of this deck are ignored and it is labeled totally as the basic 
running deck. The categories of data used are described in detail in 
previous sections (for instance, Basic Data is explained in Section 
3 .2, etc.). 
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Observation Observation 

Set Number Card Tyne Field 1 Field 2 Field 

(Column 22) (Column 77) 


CO 
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Table 10. Observation Tynes and Indicators 








































3.6.1 Trajectory Prediction Deck Arrangement 

When the function to be performed is a single trajectory pre¬ 
diction (function code ITIN = 3), the arrangement is as in Figure 5. 



TRAJECTORY PREDICTION DECK ARRANGEMENT 

When several trajectory predictions are to be done in one run 
(function code ITIN = 333, etc.) then the deck arrangement is merely 
an extension of that for the single case (Figure 6). 

3.6.2 Observation Generation Deck Arrangement 

Figure 7 and Figure 8 show the correct data setup for a single 
data generation run (function code ITIN = 4) and several cases (func¬ 
tion code ITIN = 4.4, etc.) respectively. 


74 






























Figure 6. 

MULTIPLE TRAJECTORY PREDICTIONS DECK ARRANGEMENT 


75 











































Figure 7. 

DECK ARRANGEMENT FOR 
OBSERVATION GENERATION 
SINGLE CASE. 
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Figure 8. 

Deck Arrangement for 
Two Cases of Observa¬ 
tions Generation 
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3.6.3 Orbit Determination Deck Arrangement 

There are several different deck setups that are pernrissible 
for an orbit determination run (function code ITIN = 12) according to 
the choice of several suboptions. 

Figure 9 is the simplest deck setup and applies when the ob¬ 
servation data is being input on tape (i.e., IBCDI input variable is 
used) . 

Figure 10 is the ordinary deck setup for a single case when 
the observation data is attached to the deck. Figure 11 is the deck 
setup for two such stacked cases. 

Figure 12 is the appropriate data arrangement for the situa¬ 
tion of multiple satellites (in this case three satellites). 



Deck Arrangement for Orbit Determination with Input Tape* 

Note that the BCD tape of observations must contain as its first 
set of cards the station cards followed by a TS card and then the 
observations in proper flocks followed by TF cards with a TR card 
replacing the TF card cf the last flock. 
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Figure 10. 

Deck Arrangement for Single 
Orbit Determination Run with 
All Data on Cards 
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* Note that the word "tracking ' is used interchangeably with the words 
"Orbit Determination". 
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Figure 11. 
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Figure 12. 

Deck Arrangement for Multiple 
Satellite Orbit Determination 
Run (3 Satellites) 
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3.6.4 Sample Input Decks 

Figures 13, 14, and 15 are listings of sample data decks for 
trajectory prediction, observation generation, and orbit determination. 
Besides giving the user an idea of what a typical deck looks like, 
they may also be used to familiarize the user with the operation of 
TRACE-D before he creates his own input. Another reason for running 
a test case of this sort is to determine the status of the program 
tape, for it is possible for the contents of a tape to deterioriate 
and no longer be useful 1 . 

It is not to be assumed that these sample cases exercise all 
the possible options of TRACE-D, but rather that they exhibit typical 
input n d output. A complete description of the output obtained from 
these runs is included in Section 4. 


1 In this instance, there is no cause for alarm. MITRE has the capa¬ 
bility to create a new and correct program tape in 24 hours. 
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Figure 13 (cont. 
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Figure 15 (cont 
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Figure 15 (cont. 
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Figure 15 (cont. 
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Section 4. 


SAMPLE OUTPUT 

4.0 This section presents the actual output from the test cases 

described in Section 3.6.4 whose input is listed in Figures 13, 14, 
and 15. Each test case is covered separately and the code numbers 
next to certain sections of the ouDtut correspond to the numbered 
paragraphs of the explanation that follows. 

Due to lack of space only pieces of the total output are 
shown; however, the full actual output is available. 


4.1 Trajectory Prediction Output 

The following sequenced paragraphs describe the different 
parts of trajectory prediction test case output shown in Figures 16 
to 19. 

1. The constants used in the program (c, INTEG, NUMB, 
IFLAG) are printed out here in exact input card 
image. 

2. These numbers are the geopotential inputs printed 
out exactly as read in. 

3. The last part of the input data printed out is all 
the trajectory-only input. 

4. Program name and identification of the particular 
Aerospace version. 

5. The Epoch time for this run. 

6. Indicators printed out givine the logical flow 
from segment to segment and old and new func¬ 
tion numbers. 

7. Card image of last END card in the deck. 

8. Trajectory initial conditions in three coordi¬ 
nate systems. Initial-condition values as 
shown are the result of transformations which 
have been applied to the input values. The 
transformation for the input coordinate set 
(a. 6, 6, A. r, v in this case) consists of 
conversion from decimal to octal numbers. 
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conversion of units from feet, degrees, and 
seconds to earth radii, radians, and minutes, 
and performance of corresponding inverse con¬ 
versions for outnut. The two other types of 
elements sets also require accomplishment of 
coordinate-system transformations in addition 
to the number- and units-svsterns conversions 
noted above. Accuracy of the values as 
printed therefore is subject to numerical 
roundoff errors. 

Ouantities in the left-hand column are 
position and velocity components in the basic 
vernal equinox coordinate system, with units 
of feet and of feet per second. The center 
column gives the usual ADBARV spherical sys¬ 
tem coordinates (i.e., Type-2 initial condi¬ 
tions in units of feet, feet per second, and 
degrees). The right-hand column from top to 
bottom contains orbit semi-major axis, eccen¬ 
tricity, inclination, right ascension of as¬ 
cending node, argument of perigee, and time 
of last perigee passage in minutes from mid¬ 
night of epoch day. Other units are feet and 
degrees. 

9- Identification of atmosphere model to be used 
in computing drag force. LOCKHEED and ARDC 
1959 refer to the Lockheed-Jacchia and the 
ARDC 1959 model atmospheres, respectively. 
Reciprocal of ballistic coefficient, CpA/W. 
Ballistic coefficient, W/CpA. d\ and d 2 
are values of certain constants in the 
Lockheed-Jacchia atmosphere density expres¬ 
sions. 

10. The product GM, or the universal gravita¬ 
tional constant times the mass of the earth 
(frequently designated by y), expressed in 
units of earth radii cubed per minute sauared- 
Unitless coefficients of the zonal harmonic 
terms in the earth potential field expansion 
(i.e., J 2 through J 15 ). 

Coefficients and longitudinal arguments of 
tesseral and sectorial terms in earth poten¬ 
tial field expansion, with J indicating a 
coefficient and L an argument in degrees. 

The digits following J or L are the 
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associated Legendre polynomial degree and order, 
respectively. 

11. Names of solar-svstem bodies included in computa¬ 
tion of perturbative accelerations. In this case 
none are used. 


12. Trajectory method indicator. The Cowell formula¬ 
tion of the equations of motion is the only tra¬ 
jectory method available in the present TRACE-D 
program. 


Numerical integration parameters. The Gauss- 
Jackson method (subroutine DE02) is the only 
integration method available in the present 
TRACE-D program. 

13. The indices N]_ and N 2 input in the IOB 

array respectively indicate the highest-degree 
zonal term and the highest-degree and -order 
tesseral term. No term higher than those indi¬ 
cated by these input numbers will be computed, 
regardless of the model coefficients and argu¬ 
ments that are input to the physical constants 
region. 


14. Physical constants and units conversion factors. 
Decimal and octal equivalents are given for the 
following quantities. 


a. OMEGA E(co e ) 

b. ALPHA G (ct g ) 

c. FT/KM 

d. E.R/A.U. 

e. 1-0, DISTANCE 

f. G 

g. DEG/SEC// 
RAD/MIN 


Earth rotation rate in 
radians per minute. 

Right ascension of Greenwich 
at midnight of the day of 
epoch in radians, 

Feet Per kilometer, 

Earth radii per astronom¬ 
ical unit. 

Input/output units conver¬ 
sion factor in feet per 
earth radius. 

Acceleration due to gravity 
in feet per second squared. 

Units conversion factor in 
degrees per second/radians 
per minute. 
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h. 

RELATIVE 

MASS (SUN) 

Mass of sun relative to 
mass of earth. 

i 

(VENUS) 

Mass of Venus relative to 
mass of earth, 

j. 

(JUPITER) 

Mass of Jupiter relative 
to mass of earth. 

k. 

GM 

Earth gravitational constant 
(p) in earth radii cubed per 
minute squared. 

1 . 

EARTH 

RADIUS-FT 

Number of feet per earth radius 
for internal calculations 
(rather than 1-0), 

m. 

N.M/E.R 

Number of natuical miles per 
earth radius. 

n. 

1-0, VELOCITY 

Conversion factor for input/ 
output in feet per second/ 
earth radii per minute, 

o . 

FT/SEC// 

E.R./MIN. 

Units conversion factor in 
feet per second/earth radii 
per minute, 

P- 

1/EPS 

Earth flattening (reciprocal 
of earth ellipticity), 

q- 

(MOON) 

Mass of moon relative to 
mass of earth, 

r . 

(MARS) 

Mass of Mars relative to 
mass of earth, 

s . 

(SATURN) 

Mass of Saturn relative to 
mass of earth. 

Date and time of day (Greenwich Mean Time) with which 
quantities in the print block following are to be asso 


elated. 

16. Minutes from epoch, minutes from midnight of current day, 
system time (i.e., seconds from midnight of current day), 
and current integration step size in minutes. 

17. x, y, z, r: Componenents and magnitude of the radius 
vector from geocenter to satellite in the basic coordi¬ 
nate system in units of feet. 

18. x, y, z, v: Components and magnitude of the inertial 
velocity vector with respect to the basic coordinate 
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19. 


system in units of feet per second. 

Geodetic latitude, in degrees, of the point where the 
radius vector intersects the ellipsoidal surface of the 
earth, geographic longitude measured east from Greenwich 
in degrees, altitude of the satellite above the oblate 
earth in natuical miles, and geodetic latitude of the 
subvehicle point in degrees. All latitude quantities 
are considered positive north of the equator and nega¬ 
tive south of the equator. 

20. oc, 6, 3, A: Right ascension of satellite position, de¬ 
clination of satellite position, flight path angle, and 
inertial azimuth of velocity vector in units of degrees. 

21. Revolution number, nodal period in minutes, nodal period 
decay rate in minutes per revolution, and nodal regres¬ 
sion rate in degrees per revolution. Due to the fact 
that the nodal period is calculated by simple subtrac¬ 
tion of ascending-node crossing times, this quantity 
cannot be determined until two ascending node crossings 
have been detected. Also, since the nodal-period decay 
rate is computed by differencing the nodal-period values 
at successive ascending nodes, this rate cannot be cal¬ 
culated until three ascending nodes have been crossed. 

The nodal regression rate is computed by differencing 
values of right ascension at successive ascending nodes. 

4.2 Obse r vation Generation Sample Output 

Figures 20 to 24 contain the output from the observation genera¬ 
tion test case. There are numbered descriptions corresponding to the 
number codes on the output; however, only the output that has not 
been previously described in Section 4.1 is covered here. Therefore, 
Figures 20 and 21 are included for completeness but need no further 
description (see Section 4.1 for this information). 

1. Bias inputs for data and last END card. 

2. Sigmas. In the case of a data generation run, the 
sigma table is used for the purpose of sepcifying 
the standard deviation of random noise which is to 
be applied to the generated data. Interpretation 
of the printed output shown is otherwise the same 
as described under section for Data Generation. 

3. Station locations, latitude, longitude, and 
altitude. 
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4. Data Limit Specifications. The eight columns from 
left to right contain station identification, data 
interval, minimum elevation angle, maximum eleva¬ 
tion angle, maximum range in nautical miles, and 
data-generation start and stop times in days, hours, 
and minutes from midnight of epoch day. 

The information shown indicates that the following 
specifications have been given to the program: 

a. Observations are to be generated for station 
AA at 15-second intervals whenever the com¬ 
puted local elevation angle is above zero 
degrees. 

b. No maximum elevation angle is assumed and no 
maximum range is to be considered (data gen¬ 
erated for all ranges and for all elevation 
angles greater than zero). 

c. Station AA is assumed to be active during the 
time interval from epoch until a time 24 hours 
after midnight of epoch day. 

In the case of this particular data generation run, the 
corresponding Specification-I input items for Stations 
BB through FF are identical to those for Station AA. 
However, this need not be true In general, inasmuch as 
each station is independent of the others with respect 
to these input items. 

5. Parameter indication. The message indicates that the 
range bias parameter has been selected for Station BB. 

In the case of a data generation run, radar parameters 
are selected only for the purpose of applying biases to 
the generated data, in this case a range bias on Station 
BB data. The message alluding to the correction of para¬ 
meters therefore should be ignored. 

6 . Data Types. These headings define the format of the 
table appearing beneath them. The purpose of this table 
is to define the types of data which are to be generated 
for each station. Except for the station-identification 
column header, which is printed vertically, the headings 
appear in two horizontal rows, with the top row displaying 
the symbols RANGE on the left and LONG, on the right and 
the second row the symbols SUR.R on the left and U, V on 
the right. The individual headings are interpreted in 
accordance with the following: 
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Heading 

Description or Symbol 

Units 

RANGE 

Range 

n mi 

AZMTH 

Local azimuth angle 

deg 

ELEV. 

Local elevation angle 

deg 

R. DOT 

Range rate 

ft/sec 

P. DOT 

Range-rate difference 

ft/sec 

Q. DOT 

Range-rate difference 

ft/sec 

P 

Range difference 

ft 

Q 

Range difference 

ft 

AZ. DT 

Rate of change of local azimuth 

deg/sec 

E. DOT 

Rate of change of local elevation 

deg/sec 

R. DDT 

Second time derivative of range 

ft/sec^ 

MU. VIS 

Mutual visibility 

(indicator) 

LAT 

Latitude of sub-vehicle point 

de° 

LONG 

Longitude of sub-vehicle point 

deg 

SUR. R 

Surface range, station to sub¬ 
vehicle point 

n mi 

HIGHT 

Altitude above oblate earth 

n mi 

DOPLR 

Doppler frequency shift 

cps 

LOOK 

Look angle 

deg 

VARI 

Variances 

Same as 

corresp. 

observations 

KAPPA 

Angle between radius vector 
and local vertical 

deg 

ASPCT 

Aspect angles 

deg 

ATTEN 

Signal attenuation 

db 

X, Y, Z 

x, y, z 

n mi 

T-R, D 

Topocentric right ascension 
and declination 

deg 

HR, ANG 

Topocentric hour angle 

deg 

u, v 

Horizon-scanner angles, u, v 

deg 
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The two horizontal lines associated with each station 
are the line containing the station-identification 
symbol and the line immediately below it. An X on 
the first line of a Y on the second line indicate 
that the corresponding quantity as defined by the table 
header is to be generated. In the present example, the 
indicated data types for Station AA are range, azimuth, 
elevation, range rate, height (altitude), x, v, z, 
topocentric hour angle, and horizon-sensor angles u 
and v. 

7 . Segment transfer messages. 

8 . Rise message. The time when the satellite becomes visi¬ 
ble from a particular station (at the specified minimum 
elevation angle) is obtained by interpolation and printed 
in the manner shown along with the local azimuth angle 
for the corresponding time. 

9. Time corresponding to generated data. The time shown is 

to be associated with the data quantities appearing on 
the line to the right of the time printout and on the 
line following. The time is given both in hours (0 
through 24) and minutes of the day as identified at the 
top of the output page, as well as in minutes from start 
(i.e., minutes from epoch). In the case of the oarticu- 
lar output shown, time from start corresponds to the 
time of because the epoch chosen for the run happened 

to be midnight. 

10. Maximum elevation point. The time when the elevation 
angle reaches its maximum is obtained by interpolation 
and printed along with the corresponding values for the 
elevation and azimuth angles in the manner shown. 

11. Set message. When the time at which the elevation angle 
reaches the specified minimum from above is obtained by 
interpolation, the print shown at this position occurs. 

12. Duration message. After each pass, a message is printed 
giving the time in minutes during which the elevation 
angle was above the input minimum value and the range 
was below the input maximum value. 

4.3 Orbit Dete rmin ation Sample Output 

This section contains the output from the orbit determination 
test case. Here also, only that output which has not been previously 
described is considered, although all of the output is shown in Fig¬ 
ures 25 through 31. 
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1. The station location card quantities, refraction indi¬ 
cator, sigma set, latitude, longitude, and altitude. 

2. Each line is the printout of one observation card. 
Given are the date, station ID, both the type and 
set number for these measurements, and the three 
measurements. In this case Line 1 has a Range Rate 
measurement only. 

3. This line contains a Range measurement onlv. 

4. Indicates those parameters chosen in input to be 
estimated. 

5. In this case, the snherical set of initial condi¬ 
tions, drag coefficient, GM, J 2 and J 3 are to be 
estimated. 

6 . The number of iterations after which the program 
stops. In this case, one only. 

7. A count of the parameters to be estimated. 

8 . A count of the number of observation cards read. 

9. A total obtained by assuming three observations 
per cards. 

10. A count of all the sites entered on station loca¬ 
tion cards. 

11. The number of locations that all the data has been 
packed into and then stored on scratch tape. 

12. The number of flocks (TF and TR cards) read. 

13. Segment transfer messages. 

14. Differential correction bounds. This sequence 
corresponds to the previously described sequence 
of the parameter indications. 

15. Weighting sigmas. The observation type number 
and the value of the weighting value are shown. 

16. Definition of current orbit. For each iteration 
the current values for the parameters are used 
to compute the auantities shown 

17. Node print. Each time the integrated trajectory 
crosses the equator (as determined by interpola¬ 
tion between integration steps) the date, time 
in minutes from midnight of epoch, and the rec¬ 
tangular elements (x, y, z, x, y, z,) in units 
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of earth radii and earth radii per minute are printed. 
The numerical integration interval for each satellite 
begins at epoch and ends at the time of the latest ob¬ 
servation for that satellite. A message is printed 
each time any of the events of equator crossing, start 
of thrusting, end of thrusting, or orbit adjust are 
detected during the integration interval. 

18. Station identification and system time of residuals. 
These are to be associated with all residuals which 
appear on the same line. It should be noted that 
system time is defined as seconds from midnight of 
the current day. 

19. Residuals. These are unnormalized differences be¬ 
tween the input observations (as modified by bias 
or refraction corrections) and corresponding values 
for the same observation types computed from the 
integrated trajectory position at the observation 
times. Up to six residuals for the same time are 
printed on one line. The observation type is indi¬ 
cated in parentheses immediately following the resi¬ 
dual value in each case. Note that residuals appear 
for the unweighted as well as for the weighted obser¬ 
vations . 

Identification of the foregoing observation-type in¬ 
dicators with observation descriptions or symbols 
defined elsewhere in this report is in accordance 
with the following: 


Indicator 

Description or Symbol 

Unit 

R 

Range 

ft 

A 

Azimuth 

deg 

E 

Elevation 

deg 

TR 

Topocentric right ascension 

deg 

TD 

Topocentric declination 

deg 

HA 

Topocenttic hour angle 

deg 

GR 

Geocentric right ascension 

deg 

GD 

Geocentric declination 

deg 

U 

Horizon scanner in-plane angle 

deg 

V 

Horizon scanner cross-plane angle 

deg 

H 

Altitude (h) 

ft 
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Indicator Description or Symbol 

Unit 

X 

X 

ft 

Y 

y 

ft 

Z 

z 

ft 

R 

Range 

ft 

P 

Range difference 

ft 

Q 

Range difference 

ft 

RD 

Range rate 

ft/sec 

PD 

Range-rate difference 

ft/sec 

QD 

Range-rate difference 

ft/sec 

20 . 

Observation time. The time identified 

in Item 2 in 


of seconds from midnight is eiven in alternate form, 
wherein < he day of the month and the hour and minute of 
the dry ;i:o given r; integer:'; aiv' jeconc-ii are given to 
two tlecinia I place , 

21. RMS summary o residuals. The root-mean-sauare of the 
residuals for each type of observation from each station 
is computed by the residuals editor and the result is 
printed at this location. Included are the station iden¬ 
tification, number of residuals included in the RMS (i.e., 
the total number for that station and type minus the 
number deleted by the editor on this iteration), RMS, and 
the RMS divided by input weighting sigma. The observa¬ 
tion-type indicators are the same as those used in the 
residuals print output. Units are feet, degrees and 
seconds. Interpretation of the range-rate summary for 
station MH is that 206 range rate measurements were used 
to compute the rms which equals 1.61216 feet per second. 
RMS over sigma is equal to 26.6209 feet per second. 

22. Iteration number. TRACE-D performs the tracking, or 
orbit determination function, by computing series of 
differential corrections to the parameters selected by 
the user. The iteration number is advanced each time 
the process of computing a set is repeated. The number 
may be interpreted as an indication of the number of 
times the trajectory-integration/least-squares process 
has been performed. 

23. Observation count. The number of individual observations 
included in the current least-squares computation is 
indicated. 
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24. Convergence indicator. If the weighted RMS for all 
residuals for the current iteration is less than that 
for any previous iteration, then the fitting process 
is converging and the message shown is printed. If 
the RMS obtained on the current iteration is greater 
than the smallest RMS obtained on previous iterations, 
the message 

CURRENT ITERATION IS NOT GOOD 
(RMS = O.xxxxxxxx xx) 
is printed in this position. 

25. Current solution. If the iteration is successful (i.e., 
the overall RMS has been lowered), indicated values are 
the parameter values used in the trajectory, partial 
derivative, and residuals computations which have just 
been completed. In the case of Iteration No. 1, they 
are the input values of the parameters. 

If the iteration is bad, the words 

GO BACK TO 

are printed, and the values of the parameters which so 
far have produced the lowest RMS are recovered from 
memory and printed at this point. Parameter names are 
generally self-explanatory, with the possible exceptions 
of the multiple-vehicle elements. The satellite number 
is printed before the element name for Satellites 2 
through 6. 

26. Current solution in octal digits and machine units. 

These numbers corresponding to those described in Item 
30 are given both in the octal mode and in the units 
used for internal computations. Use of these quantities 
permits bypassing units and number-system conversions 
during input and output. 

27. RMS. This is the quantity which is to be minimized in 
the tracking or orbit determination process and is the 
root-mean-square of the normalized residuals included 
in the least-squares calculations on the current 
iteration. 

28. Corrections. The result of solving the system of normal 
equations associated with the parameter which occupies 
the corresponding position in the current solution block. 
Units are feet, degrees, and seconds. 

29. Bounds. Current values of the numbers used to limit the 
size of corrections. In general, these bounds are 
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30 . 


automatically increased on a good iteration and auto¬ 
matically decreased on a bad one. 

Bounding indicator. This message will be either 


or 


HITTING BOUNDS 


NOT HITTING BOUNDS 

The first message indicates that the magnitudes of the 
corrections have been controlled by solving the system 
in such a way that the constraint implied by the bounds 
is satisfied, and the latter that the normal equations 
have been solved without applying the bounds, 

31. Next solution. Each value is the sum of the parameter 
value given in the corresponding position under ,r current 
solution” and the associated correction. These are the 
parameter values which will be used for the next solution. 

32. Next solution in octal mode and units of earth radii, 
radians, and minutes. 

33. Predicted RMS. If the fitting process is converging in 
a completely linear fashion, this will be the RMS on 
the next iteration. The comparison of this number with 
the current RMS (see Item 32) may be used to measure the 
degree to which the process has already converged. 

34. Sigma of parameters divided by sigma of the normalized 
observations. The numbers given are the square roots 
of the diagonal elements of the inverse normal matrix. 

If certain assumptions are made about the characteris¬ 
tics of the observation set, the numbers may then be 
taken as the variances on the parameter solutions. 

35. Correlation matrix, correlation coefficients for the 
parameter set. These values are computed directly 
from the covariance matrix (i.e., the inverse normal 
matrix). Rows and columns are in the same sequence 
as the Item 34 data. 

36. This is the upper'half of the normal matrix. 

37. This is the lower half of the inverse normal matrix. 
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Figure 20. Observation Generation Output listing of Input 
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Figure 32. 


General TRACE-D Flow Chart 
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Figure 33. TRACE-D Segmentation Flow Chart. 
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Section 5. 

Programming Information 
5.1 Structures and Flowcharts 

Figures 32 and 33 are both flowcharts of the TRACE-D program; 
however, each is geared to a different level of description. Figure 
32 depicts program paths determined by the function to be performed 
with only basic operations considered, while Figure 33 presents per¬ 
tinent programming details also. 

Since the TRACE-D program is extremely large (approximately 
22,000 source language cards) and can not be fully contained in core, 
the segmentation mode of operation is employed whereby the program 
becomes several coreloads. Figure 33 shows the levels of segmentation 
employed and the makeup of each and every segment. The fact that only 
the segments along a single vertical path are in core at any one point 
in time solves the storage program. 

The number appearing outside each segment block in Figure 33 
serves to identify its level and position in the structure. Those 
subroutine names underlined are the control programs, which 
among other things, transfer control to the proper segments for each 
major function. 

5•2 Tape Usage 

The TRACE-D program employs many tapes to handle input-output 
data as well as intermediate calculations. The tape requirements for 
each option consists of a basic set common to all functions (listed 
in Table 11) plus a unique option-dependent set. All of the possible 
option-dependent tapes are listed in Tables 12, 13, 14 and those that 
are mandatory are appropriately marked. Please note that no more than 
five tapes from Table 14 may be used in one run (in addition to the 
basic set, of course). 

The above description has been included for completeness; how¬ 
ever, as previously mentioned, the user may obtain a running deck 
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for any option he wishes to choose which will set up all appropriate 
trpes for him. 

Table 11. 

Basic Set of Tapes 

TAPE NUMBER NAME 

LIB1 
LIB2 

2 

3 

4 PTAPE* 

3 


FUNCTION 
Program tape //I 
Program tape #2 
System Input Tape 
System Output Tape 
Necessary Scratch Tape 
Necessary Scratch Tape 


* This is a necessary input quantity for all runs. 


TAPE NUMBER 
7 

17 


19 

15 


Table 12. 


Trajectory Prediction 
NAME 
CTAPE 

NTAPE 


DTAPE 


Tapes 

FUNCTION 

Planetary Coordinate 
Tape 

Binary trajectory 
tape (nominal) for 
trajectory differencing 
(program generated). 

Binary difference tape 
for trajectory dif¬ 
ferencing . 


A binary trajectory 
tape generated by the 
program and saved. 
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Table 13. 

Observation Generation Tapes 


TAPE NUMBER 

NAME 

FUNCTION 

11 

ETAPE 

BCD tape of observations 
program generated. 

7 

CTAPE 

Planetary tape 

10 

IFLAG(16)* 

Scratch tape 
(for trajectory) 

Necessary tape for 

this option. 



Table 14. 



Orbit Determination Tapes 

TAPE NUMBER 

NAME 

FUNCTION 

11 

IBCDI 

BCD tape of observations 

12 

IBINI 

Binary tape of compacted 
observations 

7 

CTAPE 

Planetary coordinate 
tape 

9 * 


Scratch tape (for 
normal matrix) 

5* 


Necessary scratch tape 

10* 


Scratch tape 
(trajectory) 


* Necessary tapes for this option 
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5.3 Subroutines 


This section is devoted to a complete listing of every subroutine 
in TRACE-D with a functional description and other pertinent facts 
also given. In some cases two or more routines are identical in func¬ 
tion and differ only slightly in name. This is due to the segmenta¬ 
tion rule that the same routine name cannot appear in more than one 
logical path (i.e., one complete coreload). Therefore, since some 
routines were needed in more than one path only the names were changed 
not the content. The routines are listed alphabetically in Table 15 
and the information given is self-explanatory. 
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Subroutine Listing 
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Table 15 (cont. 
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5.4 Variable Lists 


In this section are listed several variable arrays used in the 
TRACE-D program. These particular variables were chosen because they 
are the most important and critical storage areas of TRACE-D. They 
contain all of the indicators and input that dictate the logical 
program flow. 
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POSITION IN ITIN LIST OF FUNCTION BEING EXECUTED ITIN 
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OTHER THAN FIRST OR LAST COMPARISON CASES 
FOR ONE AND ONLY COMPARISON CASE 
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CNDT - PARAMETER LIST. CODES ARE IN ITRCD 
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